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ABSTRACT 


Extrusion is one of the most widely used manufacturing 
processes in industries dealing with materials such as plastics, 
polymers and pharmaceutical products In order to control the 
quality of these products and improve productivity, it is 
imperative to understand clearly the underlying mechanisms and 
determine heat transfer rates, residence time distributions, and 
the thermal and flow fields, etc within the extruder. In the 
present work, a detailed numerical study of the thermal transport 
processes in the metering section of a single screw plasticating 
extruder has been carried out The results are presented for a 
non-Newtonian material following power-law fluid behaviour The 
computations using finite-difference scheme has been performed for 
a uniform barrel temperature distribution and conducting screw A 
parametric study has also been performed. 

Typically, a single screw extruder consists of a screw 
rotating uniformly inside a barrel maintained at a fixed 
temperature. At one end there is a hopper, through which plastic 
granules are fed in. The screw extruder can be divided into three 
sections, namely, (i) the feed section which is the section near 
the hopper, (ii) the heating and conveying section where the depth 
of the screw chaoinel decreases gradually and heat is added raising 
the fluid pressure and temperature, and (iii) the metering section 
where the molten material is further heated, accompauiied by 



further increase in pressure and the material subjected to high 
shear thus enhancing the mixing 

The material in the present study is low density polyethylene 
(LDPE) which is a highly viscous non-Newtonian fluid Its 
viscosity is assumed to be following the i>ower-law model Thus, 
the viscosity of the said fluid is represented as a power law 
fvinctlon of the local strain rate and an exponential fxanction of 
temperature Therefore, the temperature and the velocity profiles 
are linked by the temperature dependence of viscosity. The 
conduction within the screw body has been considered leading to 
conjugate heat transfer between the fluid and the screw which is 
made of steel Hence, the momentum equations and the energy 
equation of the fluid and the energy equation of the screw are 
numerically solved simultaneously subject to satisfying 
compatibility requirements at the screw-fluid interface and the 
continuity equation 

As the screw channel is shallow and long, and as inertial 
forces are small in comparison to viscous (since polymers have very 
high viscosity), creeping flow approximations are made for 
conservation of momentum. The effect of viscous dissipation is 
included in the energy equation of the fluid. The coordinate 
system is fixed to the screw and thus the barrel moves in a 
directin opposite to the screw rotation. The curvature effects 
have been neglected, and therefore, the screw is treated as 
unwoxind. Therefore, the cross-section of the screw channel is 
taken as rectangular. Also, the leakage across the screw flights 



is ignored The screw is modelled as a semi-infinite body(i.e. , 
semi-infinite in the direction perpendicular to the base of the 
screw channel) in a thin layer below the screw channel base while 
the rest of the screw body behaves as adiabatic This meauis that 
the present analysis is valid for large extruders having screw 
channel height to diameter ratio much less than one. 

At this point, it is to be emphasized that the basic 
objective of this study is to develop a model which most 
accurately simulates the thermal transport processes in the 
metering section of the screw extruder To achieve this goal, a 
quasi two-dimensional model followed by a more complex fully 
two-dimensional model and most complex quasi three-dimensional 
model have been modelled The aforementioned models basically 
predict fluid velocities, pressure gradients, pressures, 
temperatures and screw body temperatures. Furthermore, the models 
compute the heat input from the barrel, local Nusselt numbers and 
bulk temp>eratures along the down-channel direction 

In the quasi two-dimensional model, the flow is taken to be 
hydrodynamically developed but thermally undeveloped In the fully 
two-dimensional model, the flow is both hydrodynamically and 
thermally undeveloped. Consequently, the solution procedure for 
the latter is more involved. For the quasi two-dimensional model, 
the iterative Newton-Raphson method requiring initial guesses for 
the pressure gradients is used to satisfy the conditions on the 
total flow rates given by the integral forms of the continuity 
equation. The two most important fall-outs of the fully 



two-dimensional model are. (a) to develop an efficient and novel 
numerical technique for computing the pressure gradients 
accurately, and (b) to show that the flow becomes locally fully 
developed at a short disteince from the inlet, — thus confirming 
the earlier use of fully developed velocity profile in the quasi 
two-dimensional model as the inlet boundary condition to be not 
physically unrealistic. 

However the main limitation of both the models is that they 
do not take into consideration the thermal convections normal to 
the screw flights as well as to the channel base and thus fail to 
predict accurately the screw surface temperature which should be 
nearly equal to barrel temperature in most part of the metering 
section The necessity of correct prediction of screw surface 
temperature auid also pressure distributions leads to the 
development of the quasi three-dimensional model which, unlike in 
the quasi and fully two-dimensional models includes additionally, 
the conservation of momentum in the direction normal to the base 
of the channel as well as thermal convections in the direction 
normal to the screw flights and to the channel base. Also, the 
contribution of thermal conduction in the direction normal to the 
screw flights is included. As expected, the inclusion of the cross 
convection effects results in accurate representation of the screw 
surface temperature and pressure distribution although the 
numerical solution procedure becomes more involved requiring the 
use of SIMPLE algorithm for the solution of flow field in the 
cross-channel planes of the channel. The results show screw 
surface temperature rising rapidly to steady values equal to or 



slightly higher than the barrel temperature in the greater part of 
the metering section This trend also tallies well with the 
reported experimental results of Marshall et al (1965) and 
Palit(1972) obtained from 63 5 mm (25 inch) extruders processing 
materials such as low density polyethylene, polypropylene and 
polyvinyl chloride These results support the use of barrel 
temperature as the thermal boundary condition at the screw surface 
by some authors in earlier quasi two-dimensional models. 

Finally, a pareimetric study based on the quasi 
three-dimensional model has also been performed to see the effect 
of variations in dimensionless volumetric flow rate (<1^). 
power-law index (n) and screw speed (N) on the pressure gradient, 
pressure, local Nusselt number and bulk temperature in the 


downchannel direction of the extruder. 
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CHAPTER 1 


INTRODUCTION 


CHAPTER 1 


1.1 Introduction 


Screw extrusion is widely used for the manufacture of 
plastics, polymers, pharmaceuticals and food products. In order to 
control the quality of these products and improve productivity, it 
is necessary to understand clearly the underlying mechanisms and 
determine heat transfer rates, residence time distributions, and 
the thermal and fluid flow fields, etc within the extruder 

Typically a single screw extruder consists of a screw 
rotating uniformly inside a barrel maintained at a fixed 
temperature (Fig. 1.1). At one end, there is a hopper through which 
plastic granules are fed in The extruder can be divided into 
three sections, namely (i) the feed section (relatively deep) 
which is the section near the hopper, (ii) the heating and 
conveying section where the depth of the screw channel decreases 
gradually and heat is added raising the fluid pressure and 
temperature, and (iii) the metering section (relatively shallow) 
accompanied by further increase in pressure and the material is 
subjected to high shear, thus enhancing the mixing While the bulk 
density of the feed solids is lower them that of the melt, the 
deep feed section enables solids to be conveyed at the same rate 
as the other section of the extruder. Otherwise the extruder may 
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be ‘starved’ , which means that certain portions of the space 
between the barrel and the screw will not be filled with polymer 

The basic principle of screw extrusion is conveying solid 
materials and melting under the combined effects of barrel heating 
and viscous heat dissipation The screw extruder, in addition to 
being an efficient pump, is also a good solids conveying device 
and an excellent melt generator 

An extruder fed by solid polymer is called a ‘Plasticating 
extruder’ . In this case, the extruder is more thaui an efficient 
pump, for it must also apply the heat energy to melt the polymer. 
There are two sources of energy for this purpose • (1) The hot 
barrel of the extruder which is equipj>ed with heaters supplies 
thermal energy to polymer by conduction, and (2) through the 
rotating shaft, mechamical energy is introduced and converted by 
viscous dissipation to thermal energy 

In addition to conveying the polymer to the die and melting 
the solids, the screw extruder also mixes and homogenizes the 
melt. Indeed, one measure of quality in the product is the 
temperature uniformity and homogeneity of the melt (i.e. the 
extrudate) leaving the extruder. The plasticating extruder can 
accomplish all these tasks despite a relatively short polymer 
residence time in the extruder, thus preventing the possibility of 
thermal degradation of the polymer. 


The solids and melt are pximped forward by the relative 



movement of screw and barrel The transporting mechanism is 
similar to a nut held in wrench with the screw rotating in the 
nut. Normally, the screw advances in the nut, but when the screw 
is physically prevented from advancing, as in the case of the 
extruder, the nut will slide in the wrench In a similar way the 
material ‘slides’ between the screw auid the barrel The screw is 
turned by a motor through a gear box, but a thrust bearing 

prevents the screw from being pushed out of the barrel by the 
pressure that develops at the die Screws, usually single 
flighted, are made of steel and often hardened for protection. 
Barrels are also made of steel, and their inner surface is either 
hardened or often coated with special wear-resistant layers The 
screw ‘floats’ in the barrel with a small clearance and the 

polymer acts as a lubricant between the tip of the screw flights 
and the barrel. 

The principal variables that can be controlled during 

operation are the screw r.p m. , the barrel temperature, and the 
pressure at the dleCe.g. with the aid of a valve). The screw 
r.p.m. depends on the size of the extruder and production rate 
requirements Typical plasticating screw extruders operate mostly 
in the range of 10-200 r.p m and, depending on their size, 
extruding upto 1600 kg/hr (higher flow rates for extruder with a 
diameter above 30 cm). The limiting factor in a given size 

extruder is quality. Higher screw speed will produce higher 
throxighput, but will usually result in deterioration of quality 
(Tadmor and Klein, 1970). 



Barrel teinp>erature setting is selected according to the 
polymer extruded The temperature has to be high enough to melt 
the polymer without causing thermal degradation. 

The main advantage of modelling plasticating extrusion 
process is that an analyst can computationally predict pressure, 
temperature, residence time distributions and degree of mixing or 
homogeneity of the extrudate and select optimum operating 
conditions without performing costly laboratory experiments. 

In the present work, a detailed numerical study of the 
thermal trauisport processes in the metering section of a single 
screw plasticating extruder has been carried out. A comparative 
study of the results based on three computer models developed by 
the present author, namely quasi two-dimensional, fully 
two-dimensional and quasi three-dimensional is presented for low 
density polyethylene (LDPE) following non-Newtonian power-law 
fluid behaviour. However, the work is general enough to include 
other power-law fluids. The computations have been carried out for 
a uniform barrel temperature distribution and conducting screw. 
The plastic granules are assumed to have been completely melted at 
the inlet of the metering section of the screw extruder and 
therefore, melting has not been considered in the present 
analysis. The models discussed here can be used for design and 
optimization of screw extruders. However, in this work, no attempt 
has been made to design such extruders 


It may be worthy of note that the models presented here for 



the metering section of the plasicatlng screw extruder can also be 
used to design a ‘melt extruder’ where the pxilymer is fed in a 
molten form In such a case, the extruder contains only one 
functional zone, namely, the metering section. Melt extruders are 
used, for example, to pelletize a polymer fed straight from the 
polymerization reactor 

1.2 Review of literature 

In recent years, experimental and nvunerical studies have been 
done to understand the underlying physical phenomena in a 
screw extruder (Tadmor and Gogos, 1979, Fenner, 1979, Tadmor and 
Klein, 1970) However, most of the studies are preliminary in 
nature and entail varying degrees of approximation. Consequently, 
mauny aspects such as effect of conducting screw and trsuisverse 
thermal convection on the thermal and flow field, warrant further 
investigation In one of the pioneering studies on the flow in a 
single screw extruder, Griffith (1962) solved the governing 
equations for the fully developed (both thermally and 
hydrodynamically) flow of an incompressible fluid through a screw 
extruder with velocity and temperature profiles essentially the 
same as those in a ch« nnel of infinite width and length. The 
effects of curvature and leadcage, across the flights were ignored. 
He acknowledged the Importance of transverse thermal convection 
i.e. convection normal to the screw flights in the recirculating 
channel flow (and the very small contribution of thermal 
conduction) in their choice of thermal boundary condition at the 
screw, by setting the temperature there equal to the barrel 
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temperature Martin’s work (1969,1970) also confirmed Griffith’s 
assumption Martin also demonstrated the potentially important 
effect of leakage flow over the screw flight, assuming that such 
flow occurs, on melt temperatures in the channel It may be noted 
that because of the intense shear in the clearance between the 
screw and the barrel, a significant amount of heat may be 
convected into the channel rather than be conducted into the 
barrel or flight tip surface although the no leakage condition 
appears to be reasonable for the purposes of computing velocity 
profiles in the channel (Martin, 1970). Marshall et al. (1965) 
and Palit (1974) reported detailed experimental results obtained 
from 63 5 mm (25 inch) extruders processing materials such as low 
density polyethylene, polypropylene and polyvinyl chloride which 
showed screw surface temperatures rising rapidly along the feed 
and compression sections to steady values slightly higher than the 
barrel temperature along most of the metering section Zaunodits 
and Pearson (1969) obtained numerical solutions for a fully 
developed, two-dimensional, isothermal and non-Newtonian flow of 
polymer melts in infinitely wide rectangular screw channels. 
Rauwendaal (1985), and Lappe and Potente (1983) have discussed the 
relationship between the volume flow rate and presure for single 
screw extruders. Lidor and Tadmor (1976) and Bigg and Middleman 
(1974) have performed theoretical analyses to determine the 
residence time distribution function aind strain distribution in 
screw extruders. Tadmor aind Gogos (1980) and Fenner (1980) have 
obtained solutions for the flow of a polymeric material in the 
various sections of an extruder. Fenner (1977) also solved the 
case of temperature profile developing along the length of the 
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screw channel Elbirll and Lindt (1984) have reported the results 
of a model in which temperature was allowed to develop along the 
screw channel In these models the screw and the barrel were 
assumed to be at the same temperature Agur and Vlachopoulus 
(1982) have studied the flow of polymeric materials, which 
included a model for flow of solids in the feed hopper, a model 
for the solid conveying zone, and a model for the melt conveying 
zone Mokhtarian and Erwin (1982) developed a mathematical model 
for mixing in a single screw extruder 

Recently. Karwe and Jaluria (1988,1990) have numerically 
guialysed by finite difference techniques the flow of a polymeric 
material (following non-Newtoniaui power law fluid behaviour) in 
the metering section of a single screw extruder. While in both the 
pap)ers the screw surface is assumed to be adiabatic, in the 
former, the barrel is at a fixed temperature and in the latter, a 
given barrel temperature distribution is specified. The velocity 
profile is assvimed to be fully developed at the inlet to the 
metering zone. Their results indicate that the temperature in the 
downstream direction has a small effect on the corresponding 
velocity field, which is determined mainly by the total volume 
flow rate Under certain circumstances, the fluid temperature near 
the exit of the extruder in the vicinity of the barrel becomes 
higher than the barrel temperature indicating transfer of heat 
frcMB fluid to the barrel. A numerical procedure is presented for 
computing the residence time distribution (RTD). Screw 
characteristics are presented in terms of bulk temperature rise 


versus screw speed. 
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The studies ©f Gopalakrishna and Jaluria (1990) and 
Gopalakrlshna et al (1992) are directed mainly at simulating 
heat and mass transfer inside the channel of a single screw 
extruder, for a power-law fluid Moisture is taken as the species 
for mass transfer , since it is of particular interest in food 
processing In this case, in addition to the strain rate and 
temperature, the viscosity of the fluid is dependent also on 
moisture As an application of this analysis to starch based food 
systems, the reactive nature of the food constituents is 
incorporated by including the rate of reaction (gelatinization) 
between moisture auid food Strong viscous dissipation is found to 
arise for typical operating conditions The effect of various 
governing parameters on pressure build-up within the extruder 
channel is determined and discussed in terms of the underlying 
physical process The moisture contours indicate that starch 
gelatinization typically occurs first at the screw root Numerical 
results are validated experimentally. 

Gupta and Kwon (1990) have investigated three dimensional 
flows of non-Newtonian viscous fluids in various geometries 
including screw channel using enriched finite element. Pressure 
flow and the temperature dependence of viscosity are not 
considered in their study. The cheuinel width to depth ratio is 
taken as 2.50. Velocity vector distribution obtained with and 
without leakage at the top right corner are shown. 


A few investigators have attemped to model the flow of food 
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stuffs in single screw extruders Bruin et al (1978) have 
discussed the flow of bio-polymers in an extruder Harper (1980) 
has discussed various aspects of food extrusion in single screw 
extruders Mohamed and Morgan (1986) have presented the results 
for average heat transfer coefficients in a single screw extruder, 
for the flow of non— Newtonian food naterials. 

1.3 Objectives of the present work 


The above literature review reveals that . 

(a) The recent computer models of the metering section of 
screw extruders are quasi two-dimensional in nature, i e the flow 
is hydrodynamically developed but thermally undeveloped at the 
entrance of the screw channel which is assumed to be infinitely 
wide. 


(b) So far no attempt has been made by earlier researchers to 
take into account in their models the conjugate heat transfer i.e. 
the coupling of conduction between the fluid and the screw, which 
is the real situation. 

(c) No work considers the effect of cross convection, i.e., 
the thermal convection normal to the screw flights The crosswise 
convection effect Cein only be included in a three dimensional 
model. Most researchers have acknowledged the contribution of this 
cross-convection by setting the screw surface temperature equal to 
the barrel temperature in their quasi two-dimensional models. 
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However, it must be noted that this is an approximate boundary 
condition, although close to reality 

(d) For the reasons mentioned in (c), no researcher could 
numerically predict the screw surface temperature distribution 

Keeping the aforementioned limitations of the previous works 
in mind, the main objective of the present study is to develop a 
model which would predict the flow and thermal characteristics of 
the polymer in the metering section of a single screw plasticating 
extruder in a most realistic manner. To acheive this goal, the 
present research is directed in the following way. 

First, a quasi two-dimensional model considering conduction 
in the screw body is developed. The velocity, pressure, 
temperature, bulk temperature, local Nusselt number at various 
sections of the extruder along the screw helix are computed The 
same results are compared with those obtained on the basis of a 
fully two-dimensional model (i.e considering uniform velocity as 
well as temperature at the entrance of the metering section) also 
developed by the present author. The idea is to see whether 
assumption of fully developed velocity profile in quasi 
two-dimensional model is adequate or not It will be shown later 
in Chapter 4 that in a fully two-dimensional model, the flow 
becomes ‘locally* fully developed at a very short distance from 
the inlet. Therefore, assummption of fully developed entry 
condition is quite Justified. Finally, a quasi three-dimensional 
model (a fully three-dimensional one would result in additional 



12 


complexities in the model without appreciable difference in the 
output, thus straining tremendously the computer time and storage 
space) which assumes fully developed velocity profile at the inlet 
and finite width of the channel is developed taking into account 
the cross convection in addition to the convection in the 
downchannel direction To check the accuracy of this model, the 
screw surface temperature distribution is compared with those 
computed on the basis of quasi two-dimensional and fully 
two-dimensional models It will be shown in Chapter 6 that the 
quasi three-dimensional model predicts the screw surface 
temperature most realistically, i.e , the nature of the 
temperature distribution matches with that observed in earlier 
experimental works. 

Furthermore, a parametric study has also been performed to 
see the effect of variations in dimensionless volumetric flow 
rate, power-law index and screw speed on the pressure gradient, 
pressure, local Nusselt number and bulk temperature in the 
downchannel direction of the extruder using the quasi 
three-dimensional model. 

All results have been obtained for low density polyethylene 
(LDPE) and steel screw combination However the models are general 
enough to include other materials following non-Newtonian 
power-law fluid behavior and other screw materials 

To sum up, the present thesis aims at (i) developing a 
rigorous computer model (starting with relatively simpler ones) 
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that Includes conduction within the screw body and thermal 
convections normal to the screw flights as well as to the screw 
root, (11) coming up with efficient and novel numerical solution 
procedures in the process and (iii) finally, checking the 
consistency of the model and understanding the complex physical 
processes within the extruder by performing a parametric study 

1.4 Major assumphons 

(a) For ease of analysis, the coordinate system is fixed to 
the screw (Fig 1.2) and thus the barrel moves in a direction 
opposite to the screw rotation 

(b) The curvature effects have been neglected, amd therefore, 
the screw is treated as unwound (Fig. 1.2). 

(c) The fluid is considered to be homogeneous following 
non-Newtonian power-law behaviour and flow and heat transfer 
steady. 

(d) As the screw chauinel is shallow and long, i e., H/W « 1, 
and as convective effects are small in comparison to viscous, 
(since polymers have high viscosity) creeping flow approximations 
are made for conservation of momentum (Schlichting, 1979) 

(e) The leakage across the screw flights have been neglected. 


(f) There is good thermal contact between the melt and the 




Ffg.1.2 Geometry of the unwound 


rectangular channel . 
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metal surfaces 

(g) All properties of the fluid except viscosity are 
constant 

(h) Viscosity of the fluid is represented as a power-law 
function of the local strain rate and an exponential function of 
temperature and is independent of pressure 



CHAPTER 2 


QUASI TWO-DIMENSIONAL 


MODEL 



CHAPTER 2 


2.1 Physical description of the quasi two-dimensional model 

The melt flow in the extruder channel has been modelled in a 
quasi two-dimensional domain The principal direction of flow is 
along the helical length of the channel Because of a cross 
component of the barrel velocity, there is a transverse 
circulation taking place also The flow is elliptic in the 
direction normal to the base of the screw root. 

With suitable approximations, the three dimensional extruder 
channel is reduced to a quasi two-dimensional domain Variation in 
the direction perpendicular to the screw flight is neglected 
assuming that channel width is much longer than channel height. 
Similarly, lubrication approximation in the direction along the 
channel obviates any variation in this direction. But, the thermal 
convection in the downstream direction (it being very strong) is 
included in the energy equation of the fluid. 

Therefore, essentially the flow becomes that of 
hydrodynamically developed but thermally undeveloped. The screw 
body is considered to be conducting, i.e., heat is conducted from 
the melt to the screw body. The screw is modelled as a 
semi-infinite body, i.e , semi-infinite in the direction 
perpendicular to the base of the screw channel. This is based on 


16 



17 


the assumption that heat will penetrate from the fluid to the 
screw body (which is very thick compared to the height of the 
screw channel in a thin layer below the screw channel base while 
the rest of the screw behaves as adiabatic Fenner (1977) in a 
review of Martin’s (1970) work regarding heat transfer between 
melt flow and a metal boundary, reported that the screw surface 
behaves as though almost thermally Insulated. The screw depth in 
present observation has been taken as ten times the height of the 
chaunnel depth The smalll curvature effect has been neglected. 

The non-dimensional channel length has been taken as 250.0 
following Karwe and Jalurla (1988,1990) At the entrance of the 
metering section, the flow is considered to be unidirectional and 
developed. Thus, at the inlet, the velocity in the direction 
normal to the screw flights (u) is considered to be zero, whereas 
the developed w-velocity profile along the chcinnel is obtained by 
solving the momentum equation for isothermal fluid subject to 
satisfying the constraints on the volumetric flow rate. 

2.2 Governing equations and boundary conditions 

The simplifications to the channel flow equations are 
introduced with the aid of the lubrication approximation. In 
essence, this involves the local replacement of the actual flow in 
the parallel gap between smooth surfaces by uniform flow between 
plane parallel surfaces. As the screw channel depth is constant in 
the downstream direction, it is reasonable to apply the 
lubrication approximation in the z-direction to velocities and 
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assume that velocities are fully developed In the dovmstream 
direction, 


u=u(x,y) 

v=v(x,y) 

w=w{x,y) 


The continuity equation becomes 


dx 



The momentum equations become: 
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(2 la) 
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(2 Ic) 


( 2 . 2 ) 


(2 3a) 


(2 3b) 


(2.3c) 


where the pressure gradient is independent of x and y . 


Owing to the presence of significant thermal convection, it 
is much less reasonable to apply the lubrication approximation in 
the z-direction to temperatures So introducing lubrication 
approxinraation in the z-direction and simultaneously retaining the 
thermal convection term, the energy equation reduces to: 
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. r ( ..*■■■ * * ®T, „ ,8^T ^ S^T, ^ T du T dy ^ 

*' S' “s * '^sy * “Si’ ■ ^ "sy * 


,9u ^ av, r aw . T aw 

T (- 5 - + ^) + yz-=- + zx-s— 

xy ay ax ay ax 


(2 4) 


The next stage of simplification frequently employed is to 
apply the lubrication approximation in the x-direction. The 
channel depth is constant in this direction, except at the flight 
For the lubrication approximation to be valid, it is assumed that 
the influence of the flight is negligible and the flow may be 
treated as though the channel were infinitely wide Symbolically, 
H « W and Eqs (2.1)-{2 4) reduce to. 


u=u(y) 

{2.5a) 

v =0 

{2 5b) 

w=w{y) 

{2 5c) 

dx 

• * ^ 

X ay 

{2 6 a) 


ay 


=0 


( 2 . 6 b) 


fl ^ 

^ * P = zy 

az z ay 


{ 2 . 6 c) 


Eq.(2.2) is automatically satisfied although continuity must 
be satisfied in terms of overall channel flow rate. The flow is 
not two dimensional in the sense that there are only two non-zero 
velocity componenets, u and w, which are functions of the third 
co-ordinate y. 
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Another very useful combination of assumptions, first 
employed by Yates (1968) is provided by taking velocity profiles 
to be locally fully developed in both the dovmstream and 
transverse directions and retaining the thermal convection in the 
dovmstream direction only The energy equation becomes 


ST 

p C vf ^ = K- 
^ p Sz f 


S^T 

f ^ T SU ^ T SW 

ay2 ^Sy y=^Sy 


(2.7) 


vfhere , 



T 

yz 



(2 7a) 

(2 7b) 


In the energy equation of the fluid, viscous dissipation 
effects are strong as polymers have high viscosity Also this 
viscosity is an exponential function of temperature and a 
power-law function of the local strain rate as is shown in 
Eq. (2.8a). 

' ’ -b(T-T ) 

H ^ JL e ° (2.8a) 

' No- 
where, 




2 , 1/2 


(2.8b) 


The screw, which is assumed to be semi-infinite in the 
direction normal to the channel base, is considered to be as 
conducting up to a thin layer beyond which the screw is considered 
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to be adiabatic Let is the depth Fig 1 2 upto which the screw 

is assumed to be conducting If H . « D , then neglecting the 

d 

small curvature effect, the portion of the screw which is 
conducting is unwound The screw then tadces the shape of a 
rectangular slab Conduction in the direction perpendicular to the 
screw base is retained whereas, conduction both in the direction 
of the downstream as well as in the transverse are neglected for 
obvious reasons. So, the energy equation in the screw becomes. 


d 

dy 



=0 


(2 9) 


where K is the conductivity of the screw material It may be 
s 

noted that the screw is made of steel. 


Except viscosity, all the properties are assumed to be 
constant . 

The boundary conditions as shown in Fig (2 1) are 


at y-0. 


dJ 
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(2.10a) 


at y = H^, 


u ® 0, w » 0, 


aT 

T_ = T , K- -5-^ 

f s f ay 


aT 

s K ^ 

s ay 


(2. 10b) 
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at y = H^+H, 

“ “ ” = \z* 10c) 

at z = 0, 

u = 0, w = w. , T- = T = T. (2 lOd) 

dev f s i 

where the subscripts f refers to fluid, s refers to screw and i 
refers to the inlet conditions 

At the inlet to the heating zone, the velocity in the 
x-direction (cross velocity) is taken as zero, following Agur and 
Vlachopoulus (1982) in which it is assumed that the solid bed is 
isothermal and the pressure varies only in the downchannel 
direction Thus the flow before the entrance is only 
uni-directional . 
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X* = x/H, y* ® y/H, z* = z/H 


u» = u / V, 


bz’ 


w* = w / V, 


bz 


e 


T - 
^b(ref ) 



P* 





» • » K(ref)- -^l) 


Pe = 


p C V, H 
^ p bz 


G = 


A V. ^ 


The dimensionless equations thus obtained are, 
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(2 12a) 

(2 12b) 

(2.12c) 

(2.12d) 

(2.12e) 

(2.12f) 

(2.12g) 

(2.12h) 


(2.13a) 
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Z-MOMENTUM: 

5p* 

dz* 



(2 13b) 


ENERGY: 

2 

Pe W* ^ ^ + G 1^*1 ”"^^6 (2 14) 

SCREW: 


ENERGY: 



(2.15) 


Similarly, the boundary conditions are also obtained in the 
non-dimensional form as: 


at y* = 0 , 


30 

ayl =0 (2.16a) 


at y* = Hyh . 


u* = 0 . w* * 0 , 


ae^ K 36 

z s s 
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(2.16b) 


at y* * 1 + H^, 


u* 


\xAb2 


V. sin^ 
b 

V. cos^ 
b 


tan^, w* * 1 , 0 ^ = 1 


(2.16c) 
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at z* *0, 

w 

U* = 0, w* = = w* . , e. = 0 = 0 (2 16d) 

V, dev f s 

bz 

The constraint on the flow in the dimensionless form are, 

u* dy* = 0 (2.17a) 

Jh^H 

w* dy* 

^»d/» 

Here, the parameter represents the dimensionless 

volumetric flow rate, generally called the throiighput emerging 
from the extruder. The parameter expresses the ratio between 
the actual flow rate along the channel and the flow rate that 
would have been achieved if all the melt moved downstream with the 
same velocity as the barrel relative to the screw. It is therefore 
a measure of volumetric efficiency. 

Thus for a given screw and temperature at the barrel, the 
parameters that govern the numerical solution are Pe, G, n, q^, 
the viscosity parameter p and K^/K^. 


= Q 

^v HV. W 
bz 


(2. 17b) 


2.4 Method of solutic 
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the inlet and boundary conditions given by Eqs (2 16a)- (2 16d) and 
constraints on the flow (Eqs (2.17a) and (2 17b)) are solved by 
means of finite difference techniques The computational domain in 
the y-z co-ordinate system is discretized to obtain the finite 
difference equations The pressure gradients in both the x- aind 
z-directions are regarded as uniform over the depth of the 
channel. The computational domain is shown in Fig. 2. 2 

2.4.1 The finite difference equations 

The governing partial differential equations can be expressed 
in the following finite difference form In the subsequent 
discussions the superscripts (•) have been dropped. The x- and 
z-momentvim equations (Eqs. (2 13a-b)) are discretized and solved in 
the channel domain. The energy equation for the channel and the 
screw are discretized. The compatibility conditions of equality of 
temperature zmd continuity of heat flux at the fluid-screw 
interface are used to obtain a discretized equation for the screw 
surface. Thus the energy equation is solved in the entire domain 
of screw body and the channel. 


2.4.2 X-MOMENTUM EQUATION 


The differential equation (Eq (2.13a)) is rewritten as: 


8x 



n-1 




^1 

8yJ 


(2.18) 
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With reference to Fig 2 2, using central difference method 
the discretized form of Eq (2.18) is obtained at grid form (i,J), 


dx 


i.J 
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7 e 
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(2.19) 


Since the pressure gradient idp/dx) is constant over all i’s 
in any j location, so Instead of double subscript (i,j) for it, 
single subscript j is used . Please note that i is the grid point 
niamber in the y-direction and J is the grid point number in the 
z-direction. Eq (2.19) may further be resolved to obtain- 
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On rearranging. 
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Fig. 2.2 Computational domain and the finite difference grid for the quasi 

two-dimensional model. 
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n-1 -^G^ 

ni scetlzatlon of 7 ^ 
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(2 22a) 
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(2.22b) 


2.4.3 Z-MOMENTUM EQUATION 


The differential equation (Eq.(2.12b)) is rewritten as: 


Sz 



(2 23) 


With reference 
the discretized form 


to Fig (2.2), using central difference method 
of Eq. (2.23) is obtained at grid point (i,j), 
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The pressure gradient dp/dz is constant over all i’s in any J 
location Therefore instead of using double subscript (i,J) for 
it, single subscript J is used Eq (2 24) may further be expanded 
as 


dz 
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. n-1 
r e 
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(w 


i+l.J 



. n-1 -^e^ 
r e 


i-1/2 




“ W ) 

1 - 1 , 


(2 25) 


Rearranging Eq (2.25), 
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(2 26) 


i+1/2, J 


n-1 -pe. 


Discetization of y e 


n-1 -P0, 


Since these two terms jr e 


n-1 -^0, 


i-1/2, J 


aind r e 


i+1/2, j 


are same for the x-momentum equations, Eqs. (2.22a-b) are used 


respectively. 
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2.4.4 Energy equation 

The energy equation is discretized separately for the 
extruder chauinel domain and the solid screw domain Using the 
compatibility conditions, a discretized equation for the 
fluid-screw interface is obtained 


CHANNEL 

The differential equation (Eq (2 12c)) is rewritten as: 

30- 3^0 , vn+1 -P0- 

Pe w ^ ^ . G [r ] e " (2 27) 

For convenience, the subscript f is dropped in the following 
discussions. In the z-direction (i.e , the downstrecim direction), 
the energy equation is parabolic. So the solution is obtained by 
marching procedure in the z-direction The discretization is 
carried out by the fully implicit marching technique 

For a fully implicit condition, one may write, 

®i.J “ 

where the current and the previous marching steps are denoted by j 
auKi J-1 respectively. 

Substituting Eq. (2.27) in Eq. (2.28) 
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r 1 {d^e • 
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Further discretization is obtained using central 
technique. 
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On rearrangement of Eq (2 30) 
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n+1 -^0 
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(2 34) 


SCREW 

The energy equation in the screw is rewritten from 
Eq. (2 14)as (subscript s is dropped)* 



(2 35) 


Using central difference technique, the Eq (2 35) is 
discretized. 


" “i.j * Vi.j ‘ 

2.4.5 Discretized equation for the interface 

To obtain the discretized energy equation at the interface 
(i.e., at the screw surface separating the fluid and the screw 
body), referring to the Fig.2 2, j fluid and 

®iface-l J screw body are expressed into Taylor series form 

dropping terms beyond second order (Carnahan et. al., 1969). 



35 


IN THE CHANNEL 


From Taylor series expansion, 
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IN THE SCREW 


From Taylor series expansion. 
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i.e. , 
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The compatibility condition at the interface is given by. 
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(2 40) 


From Eq (2 27) one can write. 


a^e 


ay 


ae 


n+i -fie 


= Pe w s Gjr e 

dz 


(2.41) 


Substituting Eqs (2 41) and Eq (2.35) in Eq (2 40) and dropping 
the subscripts f and s. 


^H-l.j ~ ^i.J 


Ay* 


Ay* 


Pe w II - G y 


n+1 -fie ) 


i.J 


®i.j ■ ®i-l,J 


Ay. 


(2.42) 


Because of solid boundary at the interface of the channel and 
the screw, w = 0. 


Rearranging the terms. 
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Ks ''i+l.J 




Ay_ _ n+1 -^e 
-o^ G r e 


i.j 


(2 43) 


AYf 

where. R * 

Please note that with reference to Fig. (2.2) i=iface, 
i-l=iface-l and i+l=iface+l. 


2.4.6 Solution scheme 

The energy Eq.(2.32) is parabolic in the z-direction, i.e., 
direction along the helix Therefore, the boundary conditions are 
necessary only at z » 0 to allow narching in the z-directin, and 
thus obtain the solution in the entire domain The bovindary 
conditions at z=0 are provided in terms of the developed velocity 
profiles at uniform temperature T=Tj^ Since u is zero at the 
inlet, obtained by solving Eq. (2 26) [modified by 
substituing 6=0] by meams of an implicit finite difference scheme. 
The iterative Newton-Rapson method is used to satisfy the 
conditions on the total flow rate given by Eq. (2 17b). This 
requires initial guess for dp/dz . The iteration is terminated 
when there is negligible change in dp/dz 

At any J location downstream of the inlet, the pressure 
gradient 3p/Sx is obtained in a way similar to that described in 
the previous paragraph for dp/dz. dp/dx is obtaind by the 
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iterative Newton-Rapson nethod subject to satisfying the the 
constraint given by the Eq. (2 17a) 

The overall solution algorithm is as follows 

(1) Gueds 9p/5x, dp/dz, u, w, 6 at downstream location z. 

1 f J 

(2) Solve Eq (2.13a) and Eq (2 13b) to obtain u, w. This 
implies solution of Eq.(2.21) and Eq (2 26) in the channel domain 

(3) Check whether Eqs (2 17a-b) are satisfied or not. If yes, 
go to Step 5. 

(4) If not, improve the guess for 9p/6x, dp/dz by iterative 
Newton-Rapson method and go back to Step 2 

(5) Solve Eq. (2 14) and Eq (2.15). This implies solution of 
Eq. (2.32), Eq (2.35) and Eq. (2.43) in the entire computational 
domain. 

(6) If max |e„^„ - e ^ 

' new guess * 

(7) Else, go to Step 2. 

(8) Print u, w, 9, 3p/3x, dp/dz 

(9) Using u, w, 0 at j , compute the same at and so 

on until z * 250 following Steps # 1-8. 
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Eq (2 21) and Eq (2 26) are solved by an Implicit finite 
difference scheme The convergence criterion used for the pressure 
gradients is 


max |^A(ap/ax), A(dp/dz) 

where A stemds for the absolute fractional change between two 
consecutive iterations. This particular convergence criterion is 
not suitable when the values of the pressure gradients become 
small, i.e , close to zero Under such circumstances, only the 
absolute change in the values of the pressure gradients is 
considered for convergence 

Eq. (2.32) is solved by unconditionally stable implicit scheme 
having an accuracy of the order (Az*) and (Ay*)^. However, to 
avoid numerical oscillations (characteristics of stable methods 
with rather large grid size in the parabolic dlection) over a very 
short distance near the inlet, since there is a step change in the 
barrel temperature at the inlet from Tj^ to Tj^, the barrel is 
assumed to be adiabatic very near the inlet 

The adiabatic boundary condition at y = 0 (Fig. 2 l)is 
represented by a three point forward difference scheme. 

The computer program is tested with a wide range of governing 
parameters. A value of 0.3 for q^ is chosen as q^=0.4 ~ 0.5 
corresponds to drag flow situation (no die situation) , where 


s 10 
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negative or no pressure gradient is generated Small values of 
represent a die of smaller opening, which restricts the flow 
considerably Therefore, q^= 0 3 is a resonable value to assume 



CHAPTER 3 


FULLY TWO-DIMENSIONAL MODEL 



CHAPTER 3 


3.1 Physical description of the fully two-dimensional model 

In the present analysis, a fully two-dimensional domain is 
considered (as shown in Fig 3.1) Variations in the transverse 
channel direction normal to the flights have been ignored because 
the width is assumed to be much larger as compared to the depth of 
the channel Following the well known lubrication approximation 
(Schlichting, 1979), the momentum equation in the y-direction can 
be omitted altogether because the component v is very small with 
respect to u and w The reason (also valid for quasi 
two-dimensional model) why lubrication approximation cam be used 
in polymer processing in which ‘films’ of several orders of 
magnitude thicker than very very thin oil films are encountered is 
the fact that viscosity of polymeric melts is also several orders 
of magnitude higher than the viscosity of the lubricating oils 

Unlike that in the quasi two-dimensional model described in 
Chapter 2, the effect of the variation in the z-direction has been 
considered Thus the model becomes fully elliptic in both the y- 

US 

and 2 -directions. The inlet w-velocity along the channel is taken 
as uniform whereas the u-velocity (i.e., the velocity component in 
the direction normal to the screw flights) is assumed to have a 
zero value following Agur and Vlachopoulous (1982). The outflow 
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boundary conditions at the exit of the channel are satisfied as 
follows For both the u- and w-velocities, the first derivative 
with respect to z has been set equal to zero because axial 
diffusion terms in the x- and z-momentum equations approach almost 
zero values at the exit of the long channel It is to be noted 
that the first and not the second derivatives are set equal to 

zero because the coefficient of viscosity, p, is not constant 

2 2 

here For the energy equation, the axial diffusion term (d T/dz ) 
at the channel outlet is set equal to zero on the assumption that 
axial diffusion is zero Also conduction within the screw body is 
considered All computations have been carried out assuming a 
steady state flow and heat transfer. 

3.2 Governing equations and boundary conditions 

The geometry of the unwound rectangular screw channel is 
shown in Fig 3 1. 

As screw channel is shallow and long, 1 e. , H/U € 1, and as 
convective effects are small in comparison to viscous, creeping 
flow approximation (Schlichting, 1979) are made for the 
conservation of momentum equations. Variations in the y- aind 
z-directions are important and considered here. The thermal 
convection is significant only along the channel length 
(z-direction) . Unlike in the quasi two-dimensional model, thermal 
conduction in the fluid and the screw are significant in both the 
transverse direction y- and downchannel direction z. It is also 
acknowledged that thermal convection normal to the screw flights 



Uniform Isothermal Flow at the Inlet 



Fig. 3.1 Boundary conditions for the fully two-dimensional model 
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(x-direction) does have considerable effect since there is a 
circulation in that direction which contributes to the mixing of 
the fluid However, the effect can be considered only in a 
three-dimensional model 

Similar to that in the quasi two-dimensional model, the screw 
body is treated as semi-infinite in the y-direction. 

With the foregoing assumptions, the momentum equations of 
the fluid and the energy equations of the fluid and the screw may 
be written as 

MELT: 

X-MOMENTUM : 


^ a 

du 'I 

^ a 

' du ] 

~ dy 

^ ay 

dz 

az] 


(3 la) 


Z-MOMENTUM : 


^ _ a_ r du 1 

dz ~ dy [ ^ dy j 



(3.1b) 


ENERGY EQUATION : 


P c 


w^ = 
p 8z 


d_ 

ay 



+ M ^ 


(3.2) 


where, *(|i) *(^) 


(3.2a) 
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In the energy equation of the fluid, viscous dissipation 
effect is strong as polymers have high viscosity This statement 
is also true for quasi two-dimensional model 




n-1 -b(T-T ) 
0 
e 


(3 3) 


where , 


y = 


t ^ \2. 

t ^ \2 

-.2 

x 2 1 

2 f — 1 + 

r + 


aw ) 

M az J 

1 sy J i 

[dzj \ 

5y J 


t1/2 


(3 3a) 


SCREW : 


ENERGY EQUATION : 



At the outlet, the fully developed flow conditions exist if 
the channel were infinitely long In practice, such a condition 
can only be approximate However, the channel length in the screw 
extruders is usually about two orders of magnitude larger than the 
other dimensions Thus it is reasonably accurate to assume fully 
developed velocity conditions at the outlet This means that the 
axial diffusion terms become smaller downstream and their 
appropach to zero are used both in the boundary conditions for the 
momentum and the energy equations at the outlet (see Eqs.(3.5b-c) 
as well as Eqs (3.11b-c)). 
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The boundary conditions as shovm in Fig 3 1 are 


at z = 0, 


u = 0, w = Q/H, = Tg = 


at z = L (exit) , 



dw 

Sz 


0 



0 



0 


at y=0. 


ai 


s 


dy 


= 0 


at y = H . ( i e. , at the interface), 
d 


u = 0, w = 0, 


ai ai 

= ^s* ^ * ‘'s ^ 


bx’ 




(3 5a) 

(3.5b) 

(3 5c) 

(3 5d) 

(3 5e) 

(3.5f) 

(3.5g) 


at y=H^ + H, 


(3.5h) 
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If the leakage across the screw flights is negligible, then 
for a total volumetric flow rate Q, the integral form of the 
continuity equation takes the form 

udy =0 (3 6a) 


r«d^« 

wdy = Q/W (3 6b} 

3.3 Non-dimensionalisation 

The governing equations are non-dimensional ised in terms of 
the following dimensionless variables, 


X* = x/H, y* = y/H, z* = z/H 


(3 7a) 


u* = u / V, , w* 
bz 


w / V, 


bz 


(3 7b) 


0 = 


T - 
^b(ref ) 



(3 7c) 



(3.7d) 



A 



(3.7e) 
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^ - *’ K(ref)- ^l) 


(3 7f) 


Pe = 


p C V, H 
p bz 


(3 7g) 


M V, 


G = 


bz 




(3.7h) 


The dimensionless equations thus obtained are 


MELT ; 


X-MOMENTUM : 


§ E * ^*1 + ^ 

5x* ■ aSy* Lr J J ^2* 



n-l 


-pe 

e 


f 

dy 


(3 8a) 


Y-MOMENTUM ; 


az* ay* 



ay J 


+ 2 


az* 



n-l 


-pe 


f ^*' 

ay _ 


(3 8b) 


ENERGY EQUATION : 


Pe w* 


ae 


f 


^ " 




+ G 



(3.9) 
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SCREW : 

ENERGY EQUATION : 

a^e a^e 

1 + 1 = 0 (3 10) 

ay*^ dz* 

Similarly, the boimdary conditions are also obtained in the 
non-dimensional form as 


at z* = 0, 


u* = 0, w* = 


Q / H _ 


W V, 


= q,, 


bz 



(3.11a) 


at z* = L / H = 250 0 (exit), 


au» 

az» 


0 


(3 11b) 


aw* 

az 


0 


(3 11c) 



(3. lid) 



(3. lie) 


at y* = 0, 



= 0 


(3. Ilf) 
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at y* = H^h 


u* = 0, w* = 0, 


(3 llg) 


ae_ K se 
f _ _s s 

^ “ K _ W 


0 = 0 - 

s f 


(3 llh) 


Compatibility Conditions 


at y» = 1 + HTH, 


V, sin^ 

“* = VAbz = »• = 1. ' 


(3 111) 


The constraints on the flow in the dimensionless form are. 


1+H^ 


u* dy* = 


(3 12a) 






•' d' 
w* dy* = 

Ju Ai 




% HV^. W 
bz 


(3 12b) 


Thus the parameters that govern the numerical solution are 


Pe, G, n, q^, the viscosity parameter p and the thermal 


conductivity ratio (K^ / K^). 
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3.4 Method of solution 

The computational domain is shown in Fig 3 2 The governing 
dimensionless Eqs (3 8)-(3 10) are solved by means of finite 
difference techniques For the sake of simplicity, in the 
subsequent discussions of this chapter, the superscript ( * ) has 
been dropped from the dimensionless equations Also, the subscript 
of e have been deleted 


3.4.1 X-MOMENTUM EQUATION 


The Eq (3 8a) is rewritten as. 


^ _ a_ 

dx ~ dy 


n-1 -pe 
r e 




n-1 -p0 
y e 



(3.13) 


n-1 -^e 

Let y e = A 


By central difference scheme, the Eq (3 13) is discretized 


as. 


5p 

A^ 

1 

- A^ 

i.l/2,J 

A ^ 

1-1/2. J 

A ^ 

1-1/2.J' 

i.J-1/2 

dx 

2 2 

Ij Ay^-" Az^ 


(3 14) 


CEnVRAL Lr?F?AR»- 

I I Ka><PUR 



Further expanding the terms. 
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ax 


J 


^1+1/2, J ^"l+l J ~ ~ ^ i-l/2,J ^^1,J ~ ^1-1, 

A 2 

AVf 


^l,J-«-l/2 ^^1 J-H 


~ ^l.J^ “ ^ i.J-1/2 


(u 


i.J 




{3 15) 


Az 


Rearranging the terms. 


ax 


■ ‘^i+l/2,j Vl.j Vl/2,j Vl..J 


{ ^i.J+1/2 j+1 ^i, j-1/2 '^i.J-1 } " 

[ ^i+1/2, j Vl/2,J ^ { ^i. j+1/2 * ^i,j-l/2 } ] "i,j 


(3 16) 


where , R = 


AYx 


Az 


3.4.2 Z-MOMENTUM EQUATION 


The Eq (3 8b) is rewritten as. 


dz 


a r 

.n-1 -pe , 

_ a r 

_n-l -pe 

^ 1 

J*ay[ 

» ' S?J 


r e 

az J 


(3.17) 


. n-1 -pe 

Let y e * A 
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By central difference scheme, the Eq (3 13) Is discretized 


as, 


dp 

dz 


J 







i-1/2 


J 


2 


dv 

^1 


i-»-l/2,J 


, 9w 


i-l/2,J 


Az 


(3 18) 


Further expanding the terms. 


az 


- ^"i-n.j ~ ^i.J^ ~ 1-1/2, J ^^i,J ~ *'l-l,J^ ^ 

« ~ A 2 

J Ay- 


^i,J-H/2 ^^i.J-H ' ^i,J^ ~ ^ i,J-l/2 ^^i,J ~ ^i,J-l^ 


Az 


(3 19) 


Rearranging the terms. 


az 


^^f ■ '^i+l/2,j ''i+l,J * Vl/2,J ''i-l.,j 


+ 2 I j+j/2 '^i, j+1 A^ j_i/2 '^i,j-l } ■ 

[ Ai+i/2.J Ai_i/2,J + 2 I Aj j^j/2 \,yi/2 } ] ”i,, 


(3.20) 



where, R = 


Az 
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3.4.3 Energy equations 

MELT : 

The energy equation of the fluid (or melt) i e , Eq. (3.9) may 
be rewritten. 


Pe w 


dd 


az 

n+1 -pe 


a^e a^e „ • 

— ~ ^ + G y e 

ay^ az^ 


(3.21) 


Let y 


= A 


Using central diffeence scheme, the Eq. (3 21) is discretized 


as. 


'■ Ay^ 


®i,>l - 201. j - 0i,j-i 
Az2 


G A 


i.J 


(3 22) 


The convective term (Pe w de/dz). , in Eq (3.22) is 
discretized using backward difference scheme because the velocity 
w in the downchauinel direction is always positive and the 
convective effects will be carried in the downstream direction. In 
other words, a first order upwind scheme is used to discretize the 
convective term, it being much more dominant than the diffusion 
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(or conduction) term in the z-direction as the Peclet number is 

large 


With this, the final discretized form is, 


”i,j V ' »ui,j - “i,j * 


2 

+ R • 


( ®i,J+l ~ ^®i.J * ®i.J-l ) 


(3.23) 


Rearranging the terms. 


®i+l.J ®i-l,J ^ ®i.J+l 




{ 


2 + 2 R + 


Pe w, , Ay, 


= - <= 


Ay^*" (3.24) 


SCREW : 

The energy equation of the screw (Eq (3 10)) is rewritten as, 




0 


(3.25) 


Using central difference scheme, the equation is discretized 
as shown in Fig. 3. 2 


6, ^ , - 20, , + 0, , . 0j ~ 20, , + 0, , , 

i+l.J 1,J 1-1, J ^ l.J-H i,J i.J-1 


= 0 (3.26) 


Ay. 


Az 
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Rearranging the terms, 


6 

-20. .*«■©. 1 

Ay w 
. ® 1 

e, , , 

- 20. . +0. . „ 

. = 0 

i+l.J 

^'’i,j '’i-i,j { 

Az J 

\ i.J+1 

^'’i,j i,j-i J 



(3.27) 


The final form becomes, 


®i+l.J 


+ e 


i-l.J 



(3.28) 


3 . 4.4 Derivation of the melt-screw interface equation 


To obtain the discretized energy equation at the interface 
(i e , at the screw surface separating the fluid and the screw 
body), referring the Fig 3.2, j fluid and 

the screw body are expanded into Taylor series form dropping terms 
beyond second order (Carnahan et al., 1969). 


CHANNEL ; 

From Taylor Series expansion. 


e 


= e 


de 


i+l,J i.J Sy 


i.J 


Ay^ + 


^1 
ay^'i.J 


Ay. 


(3.29) 


Rearranging, 
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ay 


l.J 


AVf 



i.J 


(Ay^) 

2 


(3 30) 


SCREW : 


From Taylor Series expansion. 


9 -9 - ^ 

i-i.J i.J ay 


i.J 


A + 

-^72 


Ay. 


ay 'i.J 


(3 31) 


Rearranging, 


30 

ay 


i.J 


^i.j ~ Vi,j a^e 

Ay^ 


Ay. 


ay 'i.J 


(3 32) 


The compatibility conditions obtained from Eq (3 11) are: 


ae 

ay 


f 




and. 


(3.33) 


0 


f 



ay 


and 


ae 

ay 


are substituted from Eq (3 30) and Eq.(3.32) 


respectively in Eq. (3.33) to obtain. 


®i4i,j ®i.j a^e 

^>^f ■ ay2 


Ay. 




®i-H,J ~ ®i,J 



Ay. 


2 


s 


(3.34) 
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From Eq (3 21) 



f 


= Pe w 


dz 



n+1 -pe 
-Gy e 


and from Eq (3 25) 


2 2 

a e _ d B 

2 ~ ~ 2 

ay az 


(3.35) 


(3 36) 


On substitution of Eqs (3 35) and Eq (3 36) in Eq. (3 34), 


"i.l * - ® 


%^-(Pew 


ae a^e „ • 


az2 


-Gy e 


] “Xf 

Ji.J ^ 


*^s r ®i,j ~ ®i-i, j _ a^e ^^s ] 

" ‘^f L ^^s ’az^i.j 2 J 


(3.37) 


On the solid interface, w = 0 Eq. (3 37) is further 

discretized using central difference. 


~ 01 1 I 

-iiiJ — kJ . 

Ayr 


‘ “j.J*! - 2“l,J * 
Az^ 


Gy e 


i.J J 


Ay, 


e, , - 0, . , 9, - 20, , + 0, , 4 Ay 

l.J 1-1. J _ i.J+1 iij i,J-l S 


Ay. 


Az 


(3.38) 
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On rearrangement. 


iVf I iy^ J 1-1. J 


{ 


Ay^ ^ Is '^I's 


2 (Az)^ 2 (Az)^ 


} ( 'i,>i * »i.j-i ] - 


/ '"'f . 1 . ''s 1 . h 1 . 

t (Az)2 ''I'f •'f '=f (4z)2 / ‘-J ■ 


Ay_ n+1 -p9 

G f e 


'i.J 

(3 39) 


Please note that with reference to Fig 3 2 i=iface, 
i-l=if ace-1 2 tnd i+l=iface+l. 


3.4.5 Solution scheme 

The overall solution algorithm is as follows 

(1) Guess u,w,0 at all interior grid points in the 
computational domain. Sp/5x, and 8p/3z are guessed at all 
z-locations. Note that the pressure gradients are constant over 
the height of the channel. 

(2) Solve Eq. (3.8a) and Eq. (3.8b) to obtain u and w. This 
implies solution of Eq. (3.16) and Eq. (3.20) using G-S iterative 


method using SCKl. 
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(3) With the values of u and w thus obtained from step 2, 
solve Eq (3.9) and Eq (3 10) simultaneously to obtain 0^, 0^ This 
implies solution of Eq. (3 24), Eq (3 39) and Eq (3 28) using G-S 
iterative method with SCSI to obtain 0 in the entire computational 
domain. 

(4) Check whether Eqs (3 12a) and (3 12b) are satisfied at 
each z-locatin. If yes, go to step 6. 

(5) If not, improve the guess for dp/dx, auid dp/dz Return to 
step 2 with the new values of u, w, 0 and 9p/5x, and dp/dz 

(6) Print u, w, 0-, 0 , dp/dx, and dp/dz. 

I s 

(7) Compute p, bulk temperature local Nusselt number NUjj 

(8) STOP 


The convergence criterion used for the pressure gradient is, 
max |^4(ap/ax), A(ap/az)j 5 10"* 

where A stands for the absolute fractional change between two 
consecutive iterations in the values of the pressure gradient To 
facilitate convergence, an efficient technique has been used to 
improve the guesses for the pressure gradients in subsequent 
iterations. This method is described next in the subsection 3 4.6. 



3.4.6 Method of pressure gradient correction 


Referring to Fig 3 2 , the pressure gradients 9p/9x and 9p/9z 
at the z-locations starting from the second to the last grid point 
(2, •••, N) are guessed in the beginning. First grid point (i.e , 
z = 1) is omitted because it is the inlet boundary where u and w 
velocities are known. 


Eq. (3 8a) and Eq (3.8b) are solved using the guessed pressure 
gradients In the next step, the pressure gradient 9p/9x ( ^ ) is 
perturbed by a very small amount at z = 2 Again Eq (3.8a) and 
Eq (3 8b) are solved. Now, the first column of the Jacobian 
matrix, ‘J’ shown in Eq. (3.46) is calculated The 9p/9z( i} ) is 
perturbed, Eq (3.8a) and Eq (3.8b) are solved and the second 
column of matrix, ‘J’ is obtained The same procedure is repeated 
at successive z-stations until the exit of the chaumel i.e., z = 
250 (the grid point number, N). Now the full Jacobiain matrix, ‘J’ , 
is generated. The Newton-Rapson scheme is used to get the improved 
guesses for the pressure gradients The aforementioned procedure 
is repeated in subsequent iterations till the convergence 
criterion specified in Eq. (3 40) is satisfied. The mathematical 
details are given below. 


fOp/ax,ap/«z) = udy 


wdy - 

J ur Ax 


(3.41) 


g(ap/dx,dp/dz) 


(3.42) 
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Convergence is obtained when the LHS of Eq (3 41) and 
Eq (3 42) are simultaneously zero or very close to zero 

Putting ^ = Sp/Sx, and tj = 5p/5z, Eq (3 41) and Eq. (3.42) 
respectively becomes. 


fj * *y "J ’ 


1 

udy 


i.j 


(3.43) 


«J ' "j > * 


wdy - q^ 


-■i.J 


(3.44) 


where, j = 2, 3, ••• , N 


Using Newton-Rapson scheme. 


= pg^ - [ J ] I F J- 


(3.45) 


where. 


pg*^*^= Column matrix 


•||2N-2jxl| containing new 


pressure gradients 


{( 


pg” = Column matrix ■{ |2N-2jxlj- containing old 
pressure ^ gradients 

F * Column matrix -f f2N-2] xl|- containing 


||2N-2jxl| 


functions f^, g^ to be made zero. 
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J = Jacobian matrix -I 2N - 2 


|j2N-2jx|2N-2j| 


shown in 


Eq 3 46 


"^2 

^^2 

"^2 

"^2 


arjj 


% 

% 

8g2 

% 

% 




% 

• • 




• • 



a0^ 


% 

% 

agj^ 



ai,, 




(3 46) 



CHAPTER 4 


A COMPARATIVE STUDY OF THE TWO-MODELS 



CHAPTER 4 


In the present Chapter, a comparative study of the two 

models, namely, quasi two-dimensional and fully two-dimensional 

models, developed by the present author and discussed in Chapter 2 

and Chapter 3 respectively has been carried out based on the 

results obtained by using the same for a typical non-Newtonian 

fluid (low density polyethylene) following power-law behaviour 

having n=0.5, for G=10.0, Pe=5000 0, p=l 61, q^=0.3, ^=16.54, 

Kj.=0.3 W/m-K and Kg=45 W/m-K The aforementioned input parameters 

(other than K ) have been taken from Karwe and Jaluria (1990). The 
s 

value of K has been taken from Holman (1981). 
s 

4.1 Temperature and velocity fields 

Fig 4 1 shows isotherms in the screw channel computed on the 
basis of the quasi two-dimensional euid fully two-dimensional 
model. A comparison reveals aui almost identical fluid tenperature 
distributoin in the upper half of the channel. However, in the 
lower half, the fluid temperatures for the fully two-dimensional 
case are higher than those for the quasi two-dimensional case. A 
close inspection of the isotherms in the lower half also indicates 
that in the fully two-dimensional situation, heat flows mainly in 
the negative downchannel direction which meauis that very little 
heat is conducted to the screw body. This is also confirmed by the 
isotherms in the screw body depicted iii the Fig. 4. 2 for the two 
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Fig. 4.1 Isotherms in the screw channel for the quasi 
two-dimensional and fully two-dimensional 
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H/A%A 


Fig. 4.2 Isotherms in the screw body for the quasi 
two-dimensional and fully two-dimensional 
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models Fig 4 2 also shows that the nature of the isotherms in the 
screw body for the quasi two-dimensional case is consistent with 
those in the channel (Fig 4 1). The higher fluid temperatures in 
the downchannel direction for the fully two-dimensional case may 
be attributed to the axial diffusion effect. In other words, very 
small amount of heat is conducted out into the screw body 
resulting in rise in the fluid temperature 

Fig. 4. 3 depicts the isovelocity (w*) lines in the screw 
channel for the quasi two-dimensional and fully two-dimensional 
case The isovelocity lines are almost parallel, particularly at 
the outlet. It is clear that except near the inlet, the flow 
pattern is identical for both the cases From Fig. 4. 3, it is 
obvious that for the fully two-dimensional model, the flow becomes 
‘locally’ developed at a very short distance from the inlet. This 
also shows that axial diffusion in the downchannel direction has 
very little effect on the velocity field The aforementioned 
revelations confirm that the use of fully developed inlet velocity 
profile in the quasi two-dimensional model is not physically 
unrealistic. 

Fig. 4. 4 and Fig. 4. 5 show a comparison of transverse 
temperature profiles at four downstream locations in the fluid and 
the screw body respectively for the quasi two-dimensional and 
fully two-dimensional cases. For both the models, the rate of 
change of temperature along z* is faster near the inlet and slows 
down towards the exit of the channel when ‘almost’ fully developed 
temperature profile exists. In the lower half of the channel, the 
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Fig. 4,3 Lines of constant velocity components w'*’’for 
the quasi two-dimensional and fully 
two-dimensional model. 
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2-D 

Quasi 2-D 



(q) Z*= 0.0 (b) 2*= 82.5 



Fig. 4.4 Transverse temperature profiles In the screw 
channel at four downstream locations based 
on the quasi two-dimensional and fully 
two-dimensional models. 






e 

0,0 (b) 



Fig. 4.5 Transverse temperature profiles in the screw 
body at four downstream locations based on 
the quasi two-dimensional and fully 
two-dimensional models. 
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transverse fluid temperature in the fully two-dimensional model is 
nearly viniform at the downstream locations particularly away from 
the inlet. For obvious reasons, the same trend is seen in the 
screw body Fig 4 4 also reveals that in the initial portion, the 
fluid temperature is high near the barrel, and becomes more 
uniform with increasing downstream distance In fact, as seen from 
Fig. 4 4, the temperature of the flow Increases above the barrel 
temperature, 0 > 1.0, near the exit region due to the effect of 
viscous dissipation. However, the corresponding changes in the 
velocity fields are found to be less significant, as seen in 
Fig 4.6 and Fig 4.7 

Fig. 4.6 and 4.7 depict the development of downchannel and 
cross-channel velocities, i e. , w* eind u* respectively, along the 
screw channel. While the flow becomes ‘locally’ fully developed at 
a very short distance from the inlet for the fully two-dimensional 
case, except for a slightly higher velocities in the 
two-dimensional model, there is little difference in the velocity 
profiles computed on the basis of quasi two-dimensional and fully 
two-dimensional models. However, both w* and u* velocity profiles 
in both cases keep on changing in the downstresun direction as 
because the temperature profiles are also developing in that 
direction. As in the case of transverse temperature profiles, the 
rate of development of velocity profiles is higher near the inlet 
than towards the outlet. Fig. 4. 7 also Indicates that a circulation 
exists in the x-direction as the velocity component u*’s are 
positive in the upper portion of the channel and negative in the 
lower portion. This actually creates appreciable thermal 



(c) Z* = 165.0 


(d) Z*=250 


f 

Fig. 4.6 Down channel velocity w* profiles at four 
downstream locations based on the quasi 
two-dimensional and fully two-dimensional 
models. 









(c) Z* = 165.0 


(d) Z* = 250.0 


Fig. 4.7 Crosschannet velocity u* profiles at four 
downstream locations based on the quasi 
two-dimensional and fully two-dimensional 
models. 
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convection in the direction normal to the screw flights and 
facilitates mixing of the fluid The aforesaid crosswise 
convection keeps the screw surface temperature close to the barrel 
temperature But, Fig 4 8 show that the screw surface temperature 
is lower than the barrel temperature although the fully 
two-dimensional model performs better in this respect showing 
screw surface temperature closer to barrel temperature as compared 
to the quasi two-dimensional model As discussed earlier in the 
section on Review of Literature in Chapter 1, the reported 
experimental results of Marshall et al (1965) and Palit (1974) 
indicate that screw surface temperature settles to steady values 
slightly higher than the barrel temperature along most of the 
metering section The failure of the present quasi two-dimensonal 
and the fully two-dimensonal models to predict correctly the screw 
surface temperature can be explained by the fact that the thermal 
convection in the x-direction has not been taken into account in 
the present quasi two-dimensonal and the fully two-dimensonal 
models as all variations with respect to x has been neglected 
Only a three dimensional model could consider this effect 

4.2 Pressure gradients and pressures 

Fig 4 9 and Fig. 4 10 show the corresponding pressure 
gradients ap*/9z* and the local pressure values respectively 
obtained at various downstream locations for the quasi 
two-dimensonal and the fully two-dimensonal model models. As seen 
from Fig 4 9, slightly higher pressure gradients are required in 
the case of the fully two-dimensonal model. This is because of 
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Helical Length (Z*) 


Fig. 4.8 Screw surface temperature variations along 
the helical length for the quasi two-dimen- 
sional and fully two-dimensional models. 



Pressure Gradient (0pV3z^) 
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Fig. 4.9 Variations of the dimensionless pressure 
gradients along the helical length for the 
quasi two-dimensional and two-dimensional 
models. 



Pressure (p*) 
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Fig. 4.10 Variations of the dimensionless pressure 
along the helical length for the quasi 
two-dimensional and fully two-dimensional 
models. 
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greater fluid temperature in the fully two-dimensonal case 
resulting in drop in viscosity which in turn increases the 
velocity Thus the shear rate k * increases which again increases 
the viscous drag Therefore slightly larger pressure gradients are 
required to overcome the higher viscous drag in the case of the 
fully two-dimensional model The greater negative pressure 
gradients in some part of the metering section near the inlet in 
the quasi two-dimensional case can be explained by the fact that 
although viscosity is higher because of lower temperature in the 
fluid, the shear rates are much smaller Thus the overall result 
is the larger negative pressure gradient in some portion of the 
metering section of the screw As expected, the pressure 
distribution curve for the fully two-dimensional model lies above 
that for the quasi two-dimensional case (Fig 4 10) The nature of 
the pressure curves is consistent with the pressure gradient 
distributions and tallies with Karwe and Jaluria (1990) who used 
an adiabatic screw boundary condition However, in an actual 
situation one would expect pressure to continously rise from the 
inlet of the metering section to the exit instead of a rise 
followed by a drop and then rise again This discrepancy may be 
attributed to the non-inclusion of the cross convection terms 
which can only be done in a three-dimensional model 

4.3 Bulk temperature and heat input at the barrel 

Fig 4.11 and Fig. 4. 12 show the local dimensionless bulk 
temperature and the Nusselt number NUjj at the barrel, each 
corresponding to the present quasi two-dimensional model and fully 



Nusselt Number (Nuh) 
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Fig. 4.12 Variations of Nusselt number along the 

helical length for the quasi two-dimensional 
and fully two-dimensional models. 
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two-dimensional model respectively The local bulk temperature 
^bulk Nusselt number Nu^j at the barrel are defined as 

follows 


\ulk 


w dy 


L Jh 


w dy 


®bulk 


w* 0^ dy* 


where. 




‘^in« 


f'^b(ref) ~ ^i^ 



is the heat input flux at the barrel. 


(4 1) 


(4 2) 


(4 3) 


(4 4) 


As expected, Fig. 4. 11 shows the bulk temperature in the fully 
two-dimensional model is slightly higher than that in the quasi 
two-dimensional model. In both the cases, bulk temperature exceeds 
barrel temperature near the outlet. Fig. 4. 12 shows higher Nusselt 
number near the inlet corresponding to the fully two-dimensional 
model As seen from Fig. 4. 12, the Nusselt number NUjj is very high 
near the inlet location where a steep change in the barrel 
temperature occurs. With increasing downstream distance z*. the 
value of heat flux and hence the value of Nu^, drops rapidly. 
The heat loss by the barrel is convected downstream by the flow. 
Therefore a smaller amount of heat needs to supplied to the barrel 
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downstream in order to maintain it at a constant temperature If 
the heat generation due to viscous dissipation is substantial, so 
that fluid temperature rises above the barrel temperature, then 
the heat must be removed from the barrel in order to maintain it 
at a fixed value This is indicated by the negative Mussel t 
numbers near the outlet 

4.4 The conclusions from the study of the two models 

The three most important achievements of the study of the 
quasi two-dimensional and the fully two-dimensional models are 

(i) Development of an efficient and novel numerical technique 
for computing pressure gradients accurately for the fully 
two-dimensional model 

(ii) Showing that the flow becomes ‘locally’ fully developed 
at a short distamce from the inlet, - thus confirming the earlier 
use of fully developed velocity profile in the quasi 
two-dimensional model as the inlet boundary condition to be not 
physically unrealistic 

(iil) Showing that axial diffusion has a very little effect 
on the fluid velocities and some but not significant effect on 
temperature and pressure 

However, the limitations of the models that should be kept in 


mind are : 
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(i) Fully two-dimensional model is much too computer-time 
intensive 

(ii) Both the models fail to predict realistically the screw 
surface temperature 

(iii) None of the two models shows completely realistic 
pressure distribution along the downcheinnel direction 

4.5 NEED FOR A THREE-DIMENSIONAL MODEL 

As is evident from the discussions in this chapter, the main 
limitation of both the models presented so far is that they do not 
include the thermal convections normal to the screw flights as 
well as to the channel base and therefore fail to predict 
realistically the screw surface temperature distribution which 
should be close to or slightly greater than the barrel temperature 
in most part of the metering section. It is obvious that the 
aforementioned cross-convection terms cam only be included in a 
three dimensional model. From the previous experience of the 
author with the fully two-dimensional model, it is possible to 
anticipate extremely large computer execution time and very 
complex solution procedure for a fully three-dimensional model. A 
quasi three-dimensional model on the other haind would be expected 
to give more or less similar results as those based on fully 
three-dimensional model but with much less computational effort. 
Therefore, in this study, a quasi three-dimensional model has been 
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proposed The readers are referred to the subsequent Chapters for 
a detailed discussion on the quasi three-dimensional model and a 
parametric study based on the same 



CHAPTER 5 


CXJASI THREE-DIMENSIONAL MODEL 



CHAPTER 5 


5.1 Physical description of the model 


Now that the need for a quasi three-dimensional model is 
established, the present discussion is directed towards the 
physical aspects of the said model Subsequent sections of this 
Chapter deal with the numerical modelling and the solution 
procedure. 

The primary difference between a quasi two-dimensional model 
and a quasi three-dimensional model is the inclusion in the latter 
of the cross convection terms i e , thermal convections normal to 
the screw flights as well as to the base of the screw channel. 
Obviously, in contrast with the quasi two-dimensional model, the 
channel width is finite in this case Also, thermal conduction in 
the X- and y-directions are considered Furthermore, unlike in 
two-dimewnsional analyses which neglect y-momentum following 
hydrodynamic bearing lubrication approximation, in the quasi 
three-dimensional model y-momentum is taken into account as 
transverse convection normal to the channel base arising out of 
the v-velocity is to be accoxmted for Similar to the quasi 
two-dimensional analysis, the flow is assumed to be 
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hydrodynamically developed but thermally imdeveloped That iSj the 
flow is fully developed and the temperature of the melt is uniform 
at the entrance of the metering section in the x-y plane The 
solution is obtained at the next x-y plane and so on Thus the 
solution advances in the z-direction by a suitable marching 
procedure 

Both the x-boundaries are assumed to be insulated considering 
that a negligible amount of heat is conducted between the melt and 
the uncooled screw (Fenner, 1977) Also, two-dimensional 
conduction is assumed to be present in the x-y plane of the screw 
body with both the boundaries in the x-direction acting as 
insulated considering that the screw body is sufficiently thick in 
that direction 

As before, conjugate heat transfer at the fluid-screw 
interface is modelled by assiaming the screw body (being much 
thicker compared to channel height) behaving like a semi-infinite 
solid with heat transfer taking place in a very thin layer below 
the screw surface and the rest of the screw remaining adiabatic 
As in earlier two models, the non-dimensional chamnel length is 
taken as 250 

5.2 Governing equations and the boundary conditions 


The governing dimensional equations are as follows: 
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CONTINUITY ; 


ax ay 


X-MOMENTUM 


ap _ - a f au ) ^ a r r au ^ av ) 1 
ax ■ ^ ax[ ^ax J ay [ ^ [ ay ax J J 


Y-MOMENTUM : 


ay 


„ a ( 

av ) 

^ a r 

au 

av ) ] 

^ syi 



5y 

ajJ J 


Z-MOMENTUM ; 


ap 

*^av _ 

a f 

6w 

1 . a f „ 

aw ) 

az 

ax[ 

*^ax 


dy j 

Energy equations 

of 

the 

fluid and 

the screw may 


as follows : 


FLUID : 





( ai 

ai ^ 

ax 

P c 

P 

1 ^ ax 

sy 

az 


) a 

r „ ai I 

. 9 1 

r „ ai ) 


[■^hJ 

ay 1 

i ayj 


where if> is the viscous dissipation term given by 




(5 1) 


(5 2a) 


(5.2b) 


(5.2c) 

written 


fi ^ 

(5.3) 

2 

(5.3a) 
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SCREW : 


2 2 

d T a T 

2 * 2 

8x dy 


0 


(5 4) 


The boundary conditions as shovm in Fig. 5 1 are 


at 2 = 0, 


u 


u . , V 

dev 


'^dev' 


w 


Wj . 

dev 


T 


f 



(5.5a) 


at y=0. 


ai 

dy 


s 


= 0 


(5.5b) 


at y = H , ( i e , at the interface), 
d 


u = 0, V = 0, w = 0, 


T - T r 
~ s' '"f ay 


_ „ ai 
" s ay 


(5 5c) 


at y=H^ + H, 



(5.5d) 
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at X = 0, and W 


u = 0, V = 0, w = 0, 


^ = 0 ( adiabatic condition ) C5 5e) 

The constraint on the flow in the downstream direction. 

w dy dx = Q (56) 

'*0 Jh , 

a 

5.2.1 The uncoupling of longitudinal and lateral pressure 

GRADIENTS ( PATNKAR AND SPALDING, 1972, RAITHBY AND SCHNEIDER, 

1979) 

A point to note is that the symbol p used for the pressure 
in the z-momentum Eq. (5.1c) is different from the symbol p in the 
momentum equations in x- and y-directions This is a reminder of 
the fact that in our calculation procedure, an inconsistency is 
deliberately introduced into the treatment of pressure and that 
the quantities p^^ and p are calculated separately. The pressure 
p^^ can be thought of as a form of space-averaged pressure over 
the cross-sectional pleuie. The gradient dp^^dz is calculated 
after the lateral pressure gradients are obtained. 
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5.3 Non-dimensionalisation 


The governing equations are non-dimensional ised in terms of 
the following dimensionless variables, 


X* = x/H, y* = y/H, z* = z/H 


(5 7a) 


= V / w* = w / 


(5 7b) 


T - T. 


0 = 


^b(ref) “ ^i 


(5.7c) 


V, 

p - - bz 

^ . p = ^ _ 

p 


(5 7d) 


yH- 

r* = ^. P = Mo 

bz 




A ] 


(5.7e) 


» ■ *> K(ref)- ■'l] 


(5.7f) 


Pe = 


p C V. H 
^ p bz 


(5.7g) 


G = 




(5.7h) 
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The dimensionless equations thus obtained are 


MELT : 


CONTINUITY : 


au» 

aiF ay* 


0 


X-MOMENTUM : 


ax* 



n-l au* 1 

® ax J 


a 


ay* 



'au* 



Y-MOMENTUM : 


ap* _ ^ 

ay* " ay* 



n-l ^* 1 

® ay* J 


+ 


a 


ax* 



'au* 



Z-MOMENTUM : 


ap* 

*^av 

ai^ 




^* 

ax* 


] 



(5.8) 


(5 9a) 


(5.9b) 


aw* ] 
ay* J 

(5.9c) 
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ENERGY ; 


FLUID : 


a0f 

30 

30 

* 

V* ^ 

ay* 



3^0 a^e n+1 -p0 

L + i + G y» e 

ax* ay* 


(5 10) 


SCREW : 


2 2 

3 0 3 0 

ax* ay* 


(5 11) 


The boundary conditions are also obtained in the 
non-dimensional form as: 

at 2 * = 0, 


u* = u*^ . V* = v*^ , w* = w* , , 0_ = 0^ = 0 (5.12a) 

dev dev dev f s 


at y* = 0, 


(5.12b) 



at y* = H^j/H (1 e , at the interface), 
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u* =0, V* = 0, w* = 0, 


ae. K ae 

f _ _s s 

ay* ay* 




(5 12c) 


at y* = 1 + (at the barrel surface). 


V sin^ 

“ '^hx/^hz ~ V^cosljF * V* = 0, w* * 1, 0^ = 1 (5 12d) 


at X* = 0, and V/R, 


u* =0, V* = 0, w* = 0, 


ae. 30 

= 0 { adiabatic condition ) (5 12e) 


The constraint on the flow in the downstream direction in 
the dimensionless form are. 


fH/H 

0 


1+H^ 

w* dx* dy* = 


H H V. “ 'V 
bz 


(5.13) 


5.4 Method of solution 


As mentioned earlier in section 5 2, when the pressure field 
is split in the manner indicated by Eq. (5 8)-(5.9c), two distinct 
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velocity-pressure couplings result In the first of these, 
associated with the in-plane velocity and pressure distributions, 
it is required to determine what distribution of p* should be used 
in Eq.(5.9a) and (5 9b) such that the resulting u* and v* 
velocities satisfies the equation of conservation of mass, 
Eq (5 8). For the second coupling, associated with the velocity 
distribution and pressure gradient in the parabolic direction, it 
is required to determine the value of dp*^/dz* that should be used 
in Eq (5 9c), together with the boimdary conditions, such that the 
resulting velocity distributions yields the correct dimensionless 
volumetric flow rate q^ as specified in Eq (5.13) in the 
z-direction The methods of handling these two couplings are 
presented in susequent sections of this chapter 

5AA Domain discretization 

The computational domain in the x-y plane is discretized to 
obtain the finite difference equation The pressure term does not 
appear explicitly in the momentum equations Instead, the pressure 
gradient is present and through satisfying the continuity 
equation, pressure is calculated. If the velocity components and 
pressure are at the same location and interpolated linearly for 
evaluation at the control volume faces, the resulting discretized 
equations may result in a physically unrealistic solution 
(Patankar, 1980). So to avoid such possibility, different grids 
are chosen for velocity components and pressure and they are 
suitably staggered with respect to each other • The locations where 
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the pressure and other dependent variables (temperature and 
w-velocity) are calculated, are designated as main grid points In 
practice, for computational domain discretization, the main 
control volume boxindaries are first drawn in an arbritarily 
non-uniform manner and the main grid p>oints are then placed at 
their geometric centres In this practice, the main control volume 
faces are not located midway between adjacent grid points, which 
represents less accurate finite difference representation of the 
derivatives In spite of this, the above practice is adopted here 
due to the following desirable features 

(a) The grid point at the geometric centre of the control 
volume represents the control volume better than any other point, 

(b) half control volume does not occur adjacent to a 
boundary, 

(c) discontinuities in thermophysical properties, boundary 
conditions and the source terms are easier to handle. 

The grid points corresponding to the velocity components are 
displaced with respect to the main grid points in such a way that 
they are located on the faces of the main control volume. The 
u-control volumes are staggered horizontally backwards (i e. in 
the x-direction) with respect to the main control volume, while 
the v-control volumes are staiggered vertically downwards (i.e. , in 
the y-direction). This staggering is done in such a way that the 
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displaced faces pass through the main grid points, while 
non-displaced faces lie along the main control volume faces Hence 
the resultant shape of the control volume used to compute velocity 
components and other dependent variables becomes inverted L shape 
as shown in Fig 5 2 In this figure, the main grid points are 
denoted by dots identified by letters P, W, E, S, N, SW, SE, NW, 
NE whereas u-velocity locations are denoted by arrows ( — >) 
identified by P', w', E' , s' , n', s'w', s'e', N'w' and N'e' and 


v-velocity locations are denoted by arrows ('*') identified by P" , 
W", e". S", n", s"w", S"e", n'V' and n"e" The u and v 
locations are referred to as staggered points. In the same 
Fig 5 2, the main control volume, u-control volume and v-control 
voliime are respectively identified by inclined, vertical and 
horizontal hatch lines. The type of control volume used for 
velocities near the boundaries is different as shown in Fig. 5 3. A 
complete picture of the discretized computational domain used is 
shown in Fig. 5. 4 in which NI and NJ are the total number of main 
grid points used for discretizing the computational domain, 
respectively, in the x- and y-directions In this figure, 
boundaries of the main control volume are shown as dashed lines 
The various geometrical quantities needed for discretizing the 
governing equations are also shown in Fig. 5 2, 5.3 and 5 6 in 
FORTRAN variables. Attention will now be focussed on deriving the 
discretized form of Eqs (5.8-5.11) 
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Main control 
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u-control 
Jill 1 volume 



v»control 
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Fig 5.3 Type of the control volume at the boundary. 
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Fig, 5.7 v-control volume 
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5.5 The finite difference equations 

The differential Eqs. (5 8)-(5 11) are integrated over the 
respective control volume to obtain the following finite 
difference forms. For the sake of simplicity, the superscripts (*) 
have been dropped 


5 . 5.1 X-MOMENTUM 


The x-ffiomentum equation may be rewritten from Eq (5 8a): 


9P - ? ^ 

H ^ 3x 


n-1 -80, 


au 

Sx 


] 


n-1 -80, 



(5 14) 


, n-1 -80|. 
amd let A = y e 


(5 14a) 


The Eq.(5.14) is integrated over the u-control volume shown 
in Fig. 5 6, 


rh' rP 




^ dx dy 
8x 




dx dy + 


r r» ( 

A f ^ 

^ 1 ay ax J 

j- dx dy 

(5.15) 






This may be rewritten as, 
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(p - p ) SNS(J) « 2 

P ^ 


" { ( * S ]p. ■ ( * S )„, } 

- ( * I ),, } * 

{( * H ]„,-(* S )^, } (5 16) 


The Eq. (5 16) is further discretized to obtain the final 
discretized equation. 


== ^E'^' ®u ^ ^Pw ■ Pp^ 


(5.17) 


where. 


_ 2 SNS(J) 

^E' “ SEW(I) ^uP' 


(5 17a) 


_ 2 SNS(J) * 

' SEW(IXMl) uP' 


(5 17b) 


_ 2 SEWU(I) * 

^N' DYNP(J) uP' 


(5 17c) 


®S' “ DYPSIJ) "uP' 


2 SEWU(I) 


(5 17d) 




+ a^,' + a».' 


(5.17e) 


S = A „ 
u uP 


. |v(I,J+l) - v(I-l,J+l)j- - |v(I,J+l) - v(I-l,J+l)j- j 


(5.17f) 



106 


A^' is the non-diaensional viscosity term This term is 
assumed as constant over the control volume domain and calculated 
where the variable is located. This helps in linearising the 
non-linear terms in the finite difference equation 



The Eq (5.18) is integrated over the v-control volume shown 
in Fig 5.7, 



This may be rewritten as, 
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)...) 


{(■£L-|'SU 


SNSV(J) + 


SNSV(J) 


The Eq (5.20) 


is further discre 


tized to obtain the 


discretized equation, 


3 . * 

a? Vp 


V'V' 


+ ag^^Ug^^ + 


V'^'w 

+ (Pg 


- V'V' 

- Pp) SEW(I) 


where. 


SNSV( J) . , / 

DXEP(IT vP 


SNSV(J) 4 

V' = K^TTT ^P 


V' 


2 SE W(I) . " 

DYNPVUl 


2 SEW(I) 
dypsvijT 






S = A o 
V vP 


rr , - lud.J) - ud.J 

r |u(i+i,j) - u(i+i,j' j \ 


(5 20) 

final 

(5 21) 

(5.21a) 

(5.21b) 

(5.21c) 

(5.21d) 

(5.21e) 

>}] 


✓ ^ 


(5.21f) 
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V' “ {Ps “ Pp) SEW(I) (5 2ig) 


- is the non-dimensional viscosity term which is assumed 
as constant over the control volume domain It is defined and 
calculated at the grid point of the v-velocity control volume 
This helps in linearising the non-linear terms in the finite 
difference equation. 

5.5.3 Continuity equation 


The control volume of specific interest in this case is that 
which surrounds the main grid point P in Fig 5 5 (main control 
volume in Fig. 5. 2). The continuity Eq. (5 8) is integrated between 
surfaces e, w, and n, s. The integration is performed by regarding 
those values as constant over each face to obtain. 


f” 




du , 



gjdxdy. 


'S ' 

w ^ 

S 


8v 

ay 


dx dy = 0 


(5.22) 


which yeilds, 

( Uj.- - Up- ) SNS(J) + ( ~ ^ SEW(I) * 0 (5.23) 

Before deriving the pressure-correction equation at this 
stage, let us discuss the necessity of using under-relaxation 
factor for solving the finite difference Eqs. (5 17), (5.21), 

(5.23). 
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5.6 Under-relaxation 

The generalized form of the Eq (5.17),Eq (5 21) for a control 
volume is given by, 

(aj - (5 24) 

where, ^ stands for u, v, S.. stands for S , or S ,, and 1 

(pi c c 

indicates for F>oint P , P ' respectively, for the u-momentum, and 
v-momentum. 

Eq.(5.24) appears to be linear. However, the co-efficients of 
Eq.(5.24) may themselves depend on one or more of the dependent 
variables represented by To account for the resulting 

inter-equation linkages and non-linearities, repeated solutions of 
the nominally linear form of Eq (5 24) are required Each of these 
iterative solutions is defined herein as a cycle At the 
beginning of each cycle, the co-efficients are evaluated using the 
^ values obtained in the previous cycle. With the cycle by cycle 
change in co-efficients of Eq (5.24), the resulting changes in the 
^ values can be quite large, and this may cause slow convergence 
or even divergence. To moderate the changes in successive 
solutions for and thereby improve convergence, under-relaxation 

is used. 

Patankar (1980) introduces under-relaxation into EXj.(5.24) 



through cc^ as follows, 
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Oj^ ^1 “ I ®nb ^nb ( ®i " Sj ] ^5 (5 25) 

where. Is the value of from the previous oycle, whereas 
VaiJ3oor.al and Ralthby (1984) Introduce under-relaxatlon through 
the use of an E-factor according to the following revised form of 
Eq. (5.25), 


(a, - S 


l)(l * *\*-r - =!> *°1 


(5 26) 


“l^l - I Wnb * ^ T - S,) (5 27) 


Where, = (a^ - Sj)(l + -^) 


(5.27a) 


From Eq (5 25) and (5 26), we may write. 


“r * 1 + E 


(5.28) 


Under-relaxation through the use of E-factors has better 
physical meaning than that using as discussed by VanDoormal 
suid Raithby (1984). In order to accelerate convergence, values of 
E well in excess of unity are desirable In fact values in the 
range of 1 to 20 were found to be useful for the present problems 



Ill 

5.7 Derivation of the pressure correction equation 

SIMPLEC procedure (VanDoormal and Raithby. 1984) has been 
used for handling the velocity-pressure linkages In this method 
pressure field is first guessed With this guessed pressure field, 
coefficients of the momentum equations can be evaluated allowing 
these equations to be solved to obtain the flow field. In general 
this flow field does not satisfy the continuity Eq (5.23) 
Therefore this guessed pressure field is corrected so that the 
resulting velocity field satisfies the continuity Eq (5.23) This 
is accomplised by the pressure correction equation which is 
derived in accordance with VanDoormal and Raithby (1984) by 
combining the continuity equation with truncated forms of the 
momentum equatin. After solving the pressure correction equation 
following the recommendation given in VanDoormal and Raithby 
(1984), the velocity and pressure fields are corrected, and the 
procedure is repeated until the flow field satisfies both the 
continuity and momentum equations Details of the above outline 
follow. 

Let us say that for the guessed pressure distribution p*, the 

m ^ 1 

u* and V* velocity distribution obtained , by solving the 
u-momentum and v-momentum equations respectively, satisfies 

^These asterix are different from the non-dimensional asterix. All 
the calculation are carried out with the non-dimensional values. 
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‘ I WJb * \ * V-'Pw -pp* 4- »p-V' 


(5 29) 


Dp"V*^ 


^ ^nb'^nb * ^ ~ ^jp^ * ^""^p" (5 30) 


While the u- and v-veloclties obtained from Eq (5 25) using 
the correct (but unknovm) pressure distribution p satisfy the 
continuity equation, the u* and v* velocities obtained from 
Eq.(5.29) and (5.30) do not in general satisfy the continuity 
equation Correction of the guessed pressure by p* ~ P "" P* is 
therefore necessary to correct the u* field by u' = u - u* and 
v-field by v' = v - v» The relation between p' and u' is obtained 
by subtracting Eq (5.29) from Eq (5 24) which yields to, 


' I 


hb%b 


- A^'(P' - pp 


(5.31) 


The pressure p aind velocity u that satisfy both the 
continuity and the momentum equations are. 


u = u* + u' 


(5.32) 


p = p* + p' (5.33) 

Attention will now be given to that method to find p' . The 

exact equation for p' derived from Eq. (5 31) and (5 32) and the 

continuity eqution, is complicated and unsuitable for economic 

calculation. In the SIMPLE method (Patankar, 1980), the term 

Va ,_u' in Eq. (5.31) is ignored tdiereas a term of similar 
L nb nb ^ 
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magnitude on the left hand side is retained, leading to 
inconsistency of the method 


For a ‘consistent’ approximation, leading to a suitable 
simple expression for p' , the term ®^nb^nb subtracted from 
both the sides of Eq (5 31) This yields 


( “p' - I *nb )“P' ■ ( I \b“nb - [ \b"p') * *bP' ( - Pp ] 

(5 34) 


In this method, termed ‘SIMPLEC’ by VanDoormal and Raithby 
(1984), the term | ^nb'^nb ~ E ^nb^j neglected With this 
approximation, Eq. (5 34) becomes. 


»>P' - ^nb»^ ' - V' - Pp 


l.e. , 


= dp. (p^ - p^) 


(S.36) 


where , 


d 


P 


D 



(5.37) 


Thus, 

Up^ = Up' + dp' (p^ - Pp 


(5.38) 


Eq (5 38) is an intermediate step to link velocity with 
pressure, and is called the velocity correction equation. A 
slBtilar expression can be obtained for u^'. 
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Substituting Eqs.CS.aS), (539), (541) and (542) in the 

continuity Eq (5.23), we have, 

(u*j,^ - u*pOSNS(J) + (v»j^- - v»p-)SEW(I) + 

(u^ - u^OSNS(J) + (v'- - v^-)SEW(I) = 0 (5.45) 
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This equation may be written as 


®P^P “ ® (5 46) 

where, 

\ = <1e' SNS{J) (5 473,) 

= dp- SNS(J) (5 47)j) 

^ (5 47c) 

ag = dp-- SEW(I) (5 47d) 


B = - u*p-)SNS(J) + (v»j^-- - v»p--)SEW(I)] (5 47f) 

where, dp-, d^,-, d^^". dp-- are given by Eqs (5.37), (5 40), 

(5.43), (5.44) respectively 

The direct solution on the line I in Fig. (5 4a) can be 
obtained with one application of the Thomas algorithm by suitably 
guessing the off line dependent values. Such a line by line 
solution is the basis of an iteration scheme that solves along 
each I-line and then along J-line, and repeats the pattern until 
convergence is acheived The rate of convergence of such a scheme 
depends crucially on the treatment of the off line values of the 



116 


dependent varlble. 

Let us suppose that in the current iteration, Eq(5 46) is to 
be solved along an I line, sweeping in the direction of increasing 
I According to Patankar’s suggestion (1980), the avilable 
estimate of p^ is from the previous iteratuin i e , (Pgl^. t>ut 
according to Vandoormal and Raithby (1984), the best estimate of 


Pg “ fPg^” * ‘^^^Pg " ^Pe’°^ 

Pg “ ^Pg’° 

In the expression 6 is a relaxation parameter such that for 0 

t 

= 1. p^ is taken as [p^]° This approximatation in Eq.(5.49) is 
therefore introduced into Eq (5.46) to obtain, 

(5 50) 

(5.51) 

A similar estimate Is made for solutions along J-lines. 


dpPp * ‘'eP^ * Vi * * Vi * ® 

apPj. = aj.[(p^l° ♦ (S -iHp^ - Ipp°)Hpy° 

+ Vi * Vi * Vi * ® 
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Vp ' Ve " Vw ^ Vn ^ Vs - °w 


dp. 

dz 


av 


(5 55) 


where, 

A SNS(J) 

- _ wp 

t DXEP(I) (5 55a) 

SNS(J) 

^ _ wp 

w bxpw(i) 


A SEW(I) 

- _ JfP 

N DYNP(J) 


(5.55c) 


A^ SEW(I) 
^ _ wp 

®S DYPS(J) 


(5 55d) 


* SEWd) SNS(J) 


(5 55e) 


The iinder-relaxation scheme as described in section 5.6 is 
used here to solve the Eq. (5 55). So the Eq (5 55) becomes with 
the under-relaxation method as described, 


DpWp 


ap 

E n 

^nb'^nb W dz 


1 0 
E Vp 


(5.56) 


where. Dp = ap (1 + -j^) 


(5.56a) 


Initially 


5p 


av 


5p 


is not known So, is guessed and 


dz 


subsequently it is corrected by the method described below. 
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5.8.2 PRESSURE~VEL0CITY coupling in the parabolic direction : TIC 
LONGITUDINAL PRESSURE GRADIENT 


The 


ap 


av 


dz ^ regarded as xmlform over a 

5p 

. . av 

cross-section. is obtained from the integral 

mass-conservation equation in the following manner 


At first we make an estimate of 


ap 


av 


dz 


which is denoted by 


dp* 

* 2LV 

— . This enables us to compute w» field from Eq. (5 56) The 


corresponding volume flow rate is 


• * p 

Q = 2^ Ax Ay Wj 
all P 


(5 57) 


where, all P denotes all the w-control volume in Fig 5 5 For a 

ap 

av 

correct value of — , the velocity field Wp will be found that 

yields the correct cross-sectional flow rate Q : 

Q = Ax Ay Wp (5 58) 

all P 

According to Raithby and Schneider (1979), the proposals made 

ap 

* av 

by Patankar and Spalding (1972) for finding does not 

guarantee the satisfaction of both Eq. (5.56) and Ekj. (5.58) at 
given z-locations. Though Briley’s (1974) scheme does result in an 
exact solution that requires three solutions to Eq. (5.56), Raithby 
and Scheineder (1979) found that the alternative method proposed 
by them yields the ‘exact’ result with about one-third less 
computational effort than Briley’s method. 
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2 

Motivated by the linear relation between w and — , for a 

fixed set of co-efficients, an equation for the rate of change of 


w with 


av 


8z 


is sought. Defining, 


sp, 


Sv 


R = - 


av f _ P 

dz ’ P " aR ’ 


(5 59) 


if the fp’s were known, the correct velocities would be related to 
the w’s by. 


where. 




f„AR 

(5 60) 


p p 

P 


AR = - 


- [ — 11 

(5 60a) 

[las J 

[dz J J 


value is 

chosen to 

make the total 

volume flow rate 


correct, i.e. , 


AR = 


Q - 6* 


[ Ax Ay f p 
all P 


(5 61) 


This permits both end 3 ^ to be found fro. Eq. 15.60). the 

equation for fp Is found by differentiating Eq (5.55) with respect 

^As the 00 -efflclents of the Eq.(5 56) are linearised using the 
previous Iteration values of the velocities, Eq.(5.56) beco.es 


ap. 


linear between w and 


av 
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to R to obtain, 

Apfp - Apfp + + Aj^f^ + Agfg + (5 62) 

The velocity boundary conditions have been used to obtain the 
Wp solutions so that, to maintain the same conditions on Wp, f = 0 
is applied on boundaries where the velocity is specified 

5.9 Discretization of energy equations of fluid and screw body 

5.9.1 Energy equation of the fluid 

The non-dimensional energy equation of the fluid may be 
rewritten from Eq. (5 10) as. 


Pe 




dx 


+ V 


dy 


+ w 


8z 



a^e 

dy^ 


+ G 


_ n+1 -pe 
y e 


(5.63) 


Using the continuity equation (Eq. (5 1) and 9w/5z=0 ,for 
developed velocity profile) Eq (5 63) may be rewritten in the weak 
conservation form as, 


Pe 


/a(ue) 

\ax 


a(v0) 

ay 


a(w0)\ a^0 

' “ ax2 



^ n+1 -P0 
+ G y e 


(5.64) 


Leaving only the z-direction in the left hand side, one 


obtains, 
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Pe 


a (we) 

8z 




n+1 -pe 
+ 03 - e 


- Pe 


fa(ue) 

\dx 


^ a(ve)' l 

fly I 


(5 65) 


The discretized equation is now derived by integrating 
Eq (5 65) over the main control volume shown in Fig 5 5 and over 
the interval from z and z+dz. Thus, 


Pe 


^ ^e ^+dz 

a(we) 


-z+dzji .e 


az 


dz dx dy = 


s ’’w 'z 


AX 

I 

ax[axj 

s ■'w 


dx dy dz + 


:+dz 


r 

flylflyj 


^+dzJi ,e 


dy dx dz + 


'z 'w ' s 


G V e dx dy dz 


Pe 


-z+dz 

rf 


LJ 


a(ue) 

ax 


dx dy dz - Pe 


-z+dz 



w 


a(ve) 

fly 


dy dx dz 


(5 66) 


where, the order of the integration is choosen according to the 

a (we) 

nature of the ter.. For the representation of the ter. . we 

shall assume that the grid point value of (w9) prevails throughout 
the control volu«. Fully Implicit scheiM has been used In the 
integration fro. z to z*dz. The old values at z at the grid points 
will be denoted by superscript The new values of ztdz at the 
grid points will be denoted without any superscript. Then we 

obtain from Eq, (5.66), 
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Pe[(we, - (we,“] . I^jj . 

[(^] ~ (?) ] + [ G y e ^ j SEWd) SNS(J) 

- Pe [^(ue)g - (u8)J SNS(J) - Pe [(ve)^ - (ve)J SEW(I) 

(5 67) 


The third term on the right hand side is the viscous 
dissipation term This term is assumed to be constant throughout 
the control voliime. The last two terms of Eq (5 67) are convective 
terms which are discretized following the ‘upwind scheme’ 
(Patankar, 1980) . The value of u at an Interface is equal to the 
value of u at the grid point on the upwind side of the face This 
is true for v-velocity also 


Thus, 

(ue). . [u^.o]ep - [o.-ujej. 

(ue), . - [o.-ujep 

(ve)„ = [v^.ojep - [o.-vJe„ 

(,e)^ = 


(5 68a) 


(5.68b) 


(5 68c) 


(5.68d) 


The notation Ja.bj means that the larger of a and b is to be 
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used Vhen Eq. (S 67) U repeated by this concept, the 
discretization equation becomes. 


V®P ' Ve ^ Vw * * Vs ^ 

where , 

SNs(j) r 1 

' DffiPnr ^ SNS(J) |o.-uJ 

sNS(j) r 1 

^ “ DXPW(I) * SNS(J) |uy.oJ 

% ^ DYNP(J) * SEW(I) |o,-vJ 

bYPS(j) * SEW(I) 


(5.69) 


(5 70a) 


(5.70a) 


(5 70c) 


(5 70a) 


ap + ay * a^ + as + 

.“] - [o -a 

Pe SEW(I) [[v^.o] - [o.-vj t [o.-vj - [v^.o] t 



Pe SNS(J) 


" 


Pe Wp SEW(I) SNS(J) 

— 


(5.70e) 


On rearranging, 
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Hp » ag + + Sn + ag + Pe I(u^ - u^) SNS(J) + 


(v^ - V ) SEWCI)] + 

n s 


Pe Wp SEW(I) SNS(J) 
_ 


(5 70f) 


n+1 -00 

Sp *= G f e SEWd) SNS(J) + 


Pe Wp ej SEW(I) SNS(J) 


Az 


(5.70g) 


5.9.2 Energy equation of the screw body 

The energy equation of the screw body is written as: 


2 2 

2 2 
dx dy 


(5.71) 


This equation is integrated in the main contol volume 
(Fig. 5. 5) to yield the discretized finite difference equation: 


1 



Ji 

1 

— ^ dx dy + 



Js ^ 

ax 

w 

w ^ 

s 


d^e 

dy^ 


dx dy = 0 


(5.72) 


which may be written as: 




SEW(I) = 0 (5.73) 
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On further discretization. Eq. (5 73) takes the form as 

apSp = apSp + a^e^^ + + agBg (5 74) 

where , 

, _ SNS(J) 

E " bXEPd) 

SNS(J) 

“ DXPW(I) 

_ SEW(I) 

^ ~ DYNP(J) 

„ _ SEW(I) 

^ DYNP(J) 

ap = ap + a^ + aj^ +ag (5.75e) 

5.9.3 Discretized equation at the fluid-screw interface 

In this section, a discretized equation at the fluid-screw 
interface is derived using the energy equation of the fluid and 
the screw and the compatibility conditions at the interface. The 
details of the derivation are given below (Carnahan et.al , 1969). 
The main advantage of this method is that the fluid and screw 
temperature are obtained in a single computational domain. 

Now, in the main control volume shown in the Fig 5.5 and 
Fig. (5.8), 0, , , in the fluid and 0^ . j in the screw body are 
expanded into Taylor series form dropping terms beyond the second 


(5 75a) 

(5 75b) 

(5.75c) 

(5 75d) 
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Fig. 5.8 Type of control volumes near the interace 
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order 


FLUID : 


e 


l.J+1 


0 + ^ 

i.J ay 


d^e 

2*2 
i.J ay^ 


(Ay^/z)' 


i.J 


i e. , 


^ ^i.J-*-! ~ ^i.J a^e 

8y ^ j - (Ayf/2) ay2 


Ay. 


SCREW : 


e. , . = e 


ae 


i,J-i i.J ay 


2 2 
i.J ay^ 


CAyg/2)' 


i.J 


i.e , 


80 

8y 


i.J+1 i.J + O 


i.J~"W^ 


ay 


Ay. 


i.J 


The compatibility conditions, Eq. (5 12c) are re- 
as Eq.(5.78) shovm below: 


86 

8y 


s 


_£ 5® 

K^8y 


and 0j j 



(5.76) 


(5 76a) 


(5 77) 


(5 77a) 


■written here 


(5.78) 


Substituting Eq. (5.76a) and (5.77a). one obtains 
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®i.J " 


1 

1 




ay^ 

I.J ^ J 

s 


»^f 

®i,j+i ■ ®i,j 

i _ ^ 


K 

s 


Ay^/2 

ay^ 

i.j ^ 


(5 79) 


The terms 


dy^ 


in the fluid and 


i.j 


a^e 

9y^ 


in the screw are 


ij 


substituted from Eq. (5 64) and (5 71) respectively Eq (5 64) is 
rewritten here as Eq. (5 80) and shown below 



a(ue) 

dx 


a(v0) 

Sy 


+ ^twe)] 

Sz j 


dx^ ay- 


n+1 -pe 
r e 


(5 80) 


At the screw-fluid Inteface, the velocity boundary condition 
is given by, u = v = w = 0. Therefore, the left hand side of 
Eq. (5.80) becomes zero, which yields. 




n+1 ~p0 
+ + G y e 


(5 81) 


Now energy equation of the screw is 


a^e * = 0 C5.82) 

bx dy^ 


dy^ 


(5.83) 
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Vfhen Eq (5 81) and Eq (5 83) are substituted in 
the following equation results 


®i.J " d^0 

Ay 

s 

07^ ^2 

ij ^ 


*^f 

®i,j+l ■ ®i.j . 

' a^e . -P® ' 

AVf ■ 

K 

s 

(Ay^/2) 

dx^ J 

4 


The Eq. (5 84) is integrated between e and w 
interface line: 


J [ ^ tJ 


dx = 


- 0 . 


f f r "i.j + 

s L 

tjr 


f 2 . n.l -ge ) i, 

— - + Y e I — I dx 


l^L- 


Dropping the subscripts s and f we obtain. 


®P ■ ®S _ DYPSU) 

DYFSm ® ' ' 2 


l(£l. ■ S j 


®N ■ ®P . PVWCJ) 

-mipny®™' 2 


10. • 0J 


I 


G 7 


n+1 -^0^ 

e SEW(I) 


DYNP{J) 

2 


(5 79), 


(5.84) 


over the 


(5 85) 


(5.86) 



where, 
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Ay^2 = DYPSCJ) 


Ay^/Z * DYNP(J) 


The Eq (5.86) is expressed in the form 


^P®P * ®E®E * * ®N®N * ®^S®S * ^ 


where, 

DYPS(J) *^f DYNP(J) 

®E * 2 DXEP(I) K 2 DXEP(I) 

s 


DYPS(J) ^ DYNP(J) 
* 2 DXPW(I) Kg 2 DXPW(I) 


_ *^f SEW(I) 

^ ~ K dynpTjT 

s 


_ SEW(I) 

K dyps(jT 

s 


•p “e * * ”s 


K. n-H SEW(I) DVNPCJ) 

S,^a(r . j^-^0 


(5.87) 


(5 87a) 


(5 87b) 


(5.87c) 


(5.87d) 


(5.87e) 


(5.87f) 



132 


Eq (5 87) is the desired discretized equation at the 

Interface. 

5.10 Solution algorithm 

The overall solution algorithm is as follows 

(1) Calculate u», v», w*, p* and (Sp» /dz*) for the 

av 

isothermal flow (6^=0) at the inlet in the x-y plane by solving 
Eqs. (5.8)-(5.9c) subject to satisfying Eq (5 13) 

(2) Guess u*. V*, w*, 0, (Sp* /dz*), p* in the x-y pleuie at 

3.V 

the next z-station 

(3) Solve Eq. (5.9a) and Eq. (5 9b) to obtain u* and v* in the 

said x-y plane. This implies solution of Eq. (5 17) and Eq. (5.21) 

in the x-y plane. SIMPLE algorithm of Patankar (1980) has been 
used to compute u* and v*. 

(4) Obtain the corrected pressure distribution in the plane 
using pressure correction term This means that Eq.(5.46) is 
solved to obtain the pressure correction term p' . Adding p' to p 
gives correct pressure distribution in the plane Correct u* and 
v* using Eq. (5.38) and Eq. (5.41) respectively. 

(5) Solve Eq.(5 9c) to obtain w*. This implies of solution of 


Eq. (5.55). 
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(6) Correct the pressure gradient (ap* /dz*) and w* by 

checking whether Eq.(5.13) is satisfied or not This means that 
Eq.(5 62) is solved to obtain f. Using Eq (5 60a) the pressure 
gradient /9z*) is corrected, v* is corrected by using 

Eq.(5.60). Solve Eqs. (5 69) , (5.74) , and (5 87) simultaneously to 
obtain 0 

(7) Solve Eq.{5.10) and Eq (5.11) simultiuieously to obtain 0 
This means Ek]S. (5.69), (5.87) and (5.74) are solved in the entire 
x-y domain at a particular z-station to obtain 0-dlstribution. 

(8) The maximum difference between the previous iteration 
value and the current iteration value is calculated for u*, v*, 
w*, and 0. If this value is less than a specified small value 
*eps’ , then the solution is converged Go to Step 10 

(9) Otherwise, go to Step 3. 

(10) Print u», v*, w*, 0, p*. (dp* /dz*). 


(11) Compute p*^. 0t,uik’ 

(12) If z* less than 250 0, go to Step 2 and the procedure is 
repeated for next z* locations. 


(13) STOP. 



CHAPTER 6 


RESULTS AND DISCUSSIONS 



CHAPTER 6 


In this Chapter, the results based on the quasi 
three-dimensional model developed by the present author (see 
Chapter 5) are presented The isotherms in the screw channel and 
screw body in the y-z plane, the transverse temperature (0) and u* 
and w* velocity profiles (i.e., variables as functions of y only) 
at various z* locations are compared with those obtained on the 
basis of the quasi two-dimensional model Also, velocity vector 
plots and isotherms in the x-y plane have been shown at four z* 
locations. Finally screw surface temperature, bulk temperature, 
Nusselt number, pressure gradient, and pressure along downchannel 
direction based on all three models, namely, the quasi 
two-dimensional model, the fully two-dimensional model and the 
present quasi three-dimensional model, have been compared 


It is important to note that the temperature (0), w* and u 
velocities used for comparison are the integrated average values 
over X* at each y* (see Eqs (6.1)-(6 3)) in the case of quasi 
three-dimensional model . 


0(y») 


1 

TwThT 


.W/H 

0 (x*,y*)dx* 

0 


(6 1 ) 


“ Tirar 


.W/H 

w*(x*,y*)dx» 

0 


(6 2 ) 
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u*{y") « 


fW/H) 


M/E 

u*(x*,y*)dx*‘ 

0 


(6 3) 


The input data are the same as those used for the quasi 

two-dimensional and fully two-dimensional models However, for the 

quasi three-dimensional model, the finite width (W/H) of the screw 

channel has been considered. A typical W/H equal to 5 (Palit and 

Fenner, 1972} has been used in this study. It is to be kept in 

Bind that away from the walls of the channel in the x* direction, 

the gradients of the unknowns are expected to be small Therefore, 

a non-uniform grid spacings gradually increasing in the 

x*-direction have been used. In the y*-direction, i.e , along the 

height of the channel, grid spacing is uniform The computations 

have been carried out over 23 x 17 grids and 23 x 21 grids in the 

x-y plane (i.e., the cross-section) of the chcinnel and the 

screwbody respectively. Ay_* and Ay * are 0.0667 and 0 5263 

% s 

respectively. However, a better resolution of the temperature 
gradients in the y* direction of the screw body would have been 
obtained if a finer grid spacing adjacent to the screw surface 
were taken. Az* is 2.5 as before and the computations stop at z* = 
250.0. Sufficient numerical experimentation have been performed 
for the choice of optimum number of grid points. 

6.1 Tem>erature and velocity fields 

Fig. 6.1 shows the isotherms in the screw channel for the 
quasi two-dimensional and the quasi three-dimensional model on the 
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same plot A comparison reveals that the temperature profile 
develops much faster for the quasi three-dimensional model In 
other words, the fluid temj)eratures are much greater in a larger 
portion of the channel in the downchannel direction This is 
physically expected as in the quasi three-dimensional model the 
cross convection effect has been taken into account Fig. 6 2 shows 
that the nature of the isotherms in the screw body for the quasi 
three-dimensional model is consistent with those in the channel 
Fig. 6 3 and Fig. 6. 4 show a comparison of transverse temperature 
profiles at four downstream locations in the fluid and the screw 
body respectively. The graphs in Fig. 6 3 indicate that as z» 
increases, more uniformity in the fluid temperature distribution 
is achieved in t * case of the quasi three-dimensional model This 
better uniformity can be attributed to better mixing due to the 
cross convection effect. The quasi two-dimensional model does not 
show this level of uniformity in the temperature profiles because 
of non-inclusion of the thermal convection normal to the screw 
flights in the same. Also, the fluid temperature is higher for the 
quasi three-dimensional model. However, the difference with 
respect to the quasi two-dimensional model is higher in the middle 
8ukl gradually decreases towards the exit of the channel As 
expected. Fig. 6. 4 shows that the screw temperature profiles are 
almost uniform at each z* location and the screwbody temperatures 
are higher in the case of the quasi three-dimensional model. 


Fig. 6. 5 depicts the Isovelocity (w*) lines in the screw 
channel for the quasi two-dimensional and the quasi 
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Fig. 6.2 Isotherms in the screw body for the quasi 

two-dimensional and quasi three-dimensional 
models. 
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Quasi 3-D 



Fig. 6.4 Transverse temperature profiles in the screw 
body at four downstream locations based on 
the quasi two-dimensional and quasi 
three-dimensional models. 
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three-di»ensional ■odel on the same plot A comparison reveals 
that w* velocities based on the quasi three-dimensional model are 
higher in the upper half of the channel in the middle of the 
channel length but they coincide with the quasi two-dimensional 
^•-velocities towards the channel exit. This is clearly revealed 
in Fig 6.6 This is because higher fluid temperature result in 
drop in viscosity which in turn increases the flow velocity Near 
the entry and the exit, the temperature difference is small and 
hence the w* velocities for the two models almost coincide. 

Fig. 6. 7 shows the development of u* velocity profiles along 
the downchannel direction. Vfhile for the two models, u*-velocities 
differ in the middle of the channel length as well as at the 
inlet, the former coincide at the exit The higher u*-velocities 
for the quasi three-dimensional model in the middle of the channel 
length is attributed to the higher fluid temperature resulting in 
drop in viscosity whereas near the exit the temperature difference 
corresponding to the two models Is small and hence the u--»eloclty 
profiles for the quasi two-dimensional and the quasi 
three-dimensional models become Identical. As for the deviation at 
the inlet, the reason is as follows. The inlet velocity u» has 
been taken as aero for the quasi two-dimensional case whereas for 
the quasi three-dimensional model, the u» and V and w* profiles 
are obtained by solving g,s. (5.8)-(5.9c) subject to satisfying 

*0 in Eqs.(5.9a)-(5.9c). 


Eq. (5.13) by putting 
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Fig, 6.6 Downchannel velocity w* profiles at four 
downstream locations based on the quasi 
two-dimensional and quasi three-dimensional 
models. 
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Fig. 6.7 


.sschannel velocity u* profiles at four 
vnstream locations based on quasi 
dimensional and quasi three-dimensional 





6.2 Velocity vector plots and isotherms in x-y plane of the 

CHANNEL 

Fig 6 8 shows the velocity vector plots in the x-y plane 
at four downs, t ream locations The plots clearly show the cir- 

culatory nature of the flow in the cross-sectional plane of 
the extruder This circulation results in the cross thermal 
convection enhancing the overall mixing process It is also 
noticed that v* velocities are negligible except near the right 
and left boundary walls meaning a very weak flow and thermal 
convection normal to the channel base Fig 6 9 shows isotherm 
plots in the x-y planes at four locations in the downchannel 
direction. The figure depicts visually the increase of fluid 
temperatures as the fluid, moves towards the exit of the channel 
At 2* == 12.5, i.e., near the inlet of the channel, it is seen a 
large core of relatively cold fluid sandwiched between hot fluids 
in the upper and lower region of the channel cross-section Near 
the right and left boundaries, temperatures are large because the 
sides are assumed to be insulated At higher z*, i e , at z* = 
82 5, a dramatic rise in the fluid temperature is noticed This is 
due to the fact that more heating of the fluid is now taking place 
as the fluid is exposed to larger surface of the heated barrel and 
the effect of viscous heat generation is now more pronounced The 
temperature gradient in the y-direction as well as in the 
x-direction decreases which means that compared to that at z* 

12 5. fluid is better mixed. At z*=165 0, more uniformity in the 
temperature distribution is achieved and a major portion of the 
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Fig. 6.8 Velocity vector plots in the cross-sectional 
(x-y) planes at four down stream locations 
for the quasi three-dimensional model. 










fluid (except near the screw surface) Is at temperatures higher 
than that at the barrel This is due to greater viscous 
dissipation At z**250.0. i e., at the exit of the channel, it is 
clearly seen that the entire fluid in the x-y plane is at a 
temperature higher than that at the barrel The fluid is now at 
almost uniform temperature, i e , thoroughly mixed To obtain a 
well -mixed fluid at the die is one of the main objectives Also, 
since near the exit, barrel is being heated by the fluid, 
appropiate amount of heat must be removed from the barrel to 
maintain it at a fixed temperature Another important and 
interesting observation is that at the exit the isotherms take the 
shape of streamlines (as is evident from the vector plots of 
Fig 6.8). This phenomenon revealed by the present numerical study 
matches also with the expectation of Griffith (1962) 

6.3 COMPARtSOW OF RESULTS BASED ON THE QUASI TWO-DIMENSIONAL, 

FULLY TWO-DIMENSIONAL AND QUASI THREE-DIMENSIONAL MODELS 

Finally, screw surface temperature, pressure gradient, 
pressure, Mussel t number and bulk temperature along downchannel 
direction z* based on the cjuasi two-dimensional, fully 
two-dimensional and quasi three-dimensional models have been 
compared in order to see how the cjuasi three-dimensional mcxlel 
fares In contrast with the other two models. 

Fig. 6. 10 shows that the screw surface temperature is close to 
barrel temperature {e^^-l.O) beyond z«=70.0 and greater than barrel 
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r«g. 6.10 Screw surface temperature variations along 
the helical length for all three models, 
-quasi two-dimensional, fully two-dimen- 
sional and quasi three-dimensional. 
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temperature beyond z**190.0. In other words, the screw surface 
temperature is almost equal to or greater than the barrel 
temperature in most part of the metering section This trend is 
also confirmed by earlier experimental works of Marshall et al. 
(1965) and Palit (1972) obtained from 63.5 mm (2.5 inch) extruders 
processing materials such as LDPE, polypropylene and polyvinyl 
chloride. Griffith also assiuned the screw surface temperature to 
be equal to the barrel temperature on the ground that very large 
Peclet nvuBber associated with melt flow in a screw channel implies 
that temperature is constant along a streamline If there is no 
leakage over the flight tip, the screw and the barrel boundaries 
form a single continuous streamline. 

It is clear that quasi two-dimensional and the fully 
two-dimensional models have failed to predict physically realistic 
screw surface temperature behaviour 

Fig. 6. 11 reveals that in sharp contrast with the quasi 
two-dimensional and fully two-dimensional models. pressure 
gradient computed by the quasi three-dimensional model never goes 
to negative value and therefore, pressure continually increases to 
Its final value at the exit of the channel (Fig. 6. 12) which Is 
according to the physical expectations Another re,»rlcable feature 
is that pressure predicted by the quasi three-dl«nslonal «»lel Is 
.uch greater than that by the other two aodels. Therefore, quasi 
three-dlnenslonal «del In this regard also perferas superior to 
^al twcdl^nslonal «xlel and fully two-dl»=nslonal «del. It Is 
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Fig. 6.11 Pressure gradient distribution along the 
helical length for all three models, -quasi 
two-dimensional, fully two-dimensional and 
quasi three-dimensional. 
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Fig. 6.12 Pressure distribution along the helical length 
for all three models,- quasi two-dimensional, 
fully two-dimensional ond quasi 
three-dimensional . 
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to be noted that pressure predicted by the quasi three-dimensional 
model is the average pressure (p^^) computed from Eq. (5 92) 
Fig 6 13 depicts that Nusselt number changes its sign from 
positive to negative (i.e.. the fluid temperature near the barrel 
becomes higher than the barrel temperature) at around z*=150 0 for 
all the three models Furthermore, the quasi three-dimensional 
model predicts greater Nusselt number It is to be noted at this 
point that for the quasi two-dimensional and fully two-dimensional 
models, in order to avoid numerical oscillations, the barrel 
boundary is taken as adiabatic near the inlet (till z*=10) That 
is why the Nusselt number is zero at the inlet (i.e , at z»=0) 
For obvious reason, this is not done in respect of the quasi 
three-dimensional model 

Fig 6.14 shows that the bulk temperature predicted by the 
quasi three-dimensional model is greater in most part of the 
metering channel than that predicted by the other two models 
However, at the exit of the channel, bulk temperature based on the 
quasi three-dimensional model coincides with that computed by the 
fully two-dimensional model. This may be due to the effect of 
axial diffusion considered in the fully two-dimensional 


model. 
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Fia 6 14 Bulk temperature distribution along the 
* helical length for all three models, -quasi 
two-dimensional, fully two-dimensional and 

quasi three-dimensional model. 
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has been shown in Chapter 6 that the quasi 
three-dimensional model simulates the physical processes in the 
metering section of a single-screw plasticating extruder in a most 
realistic manner Therefore, the quasi three-dimensional model can 
be used for the simulation of the performance of different size 
extruders with various screws and operating conditions to 
determine their relative performance in the desired application 
This procedure will yield an optimized design and valuable 
^^ormation on the performance of the extruder in a variety of 
operating conditions However. since the aim of this work is not to 
arrive at an optimized design, the parametric study has been 
performed on a single screw extruder thus excluding variation of 
W/H ratio. 

The purposes of the present parametric study are two fold: 
(i) to see whether the quasi three-dimensional model is 
consistently behaving in a physically realistic manner when 
operating variables are changed at a time and (ii) to gain a great 
deal of insight auid information on the process 

The effect of three variables, namely, non-dimensional 
volumetric flow rate power law index (n) and screw speed 

(N), on pressure gradient, pressure, Nusselt number and bulk 
temperature of a screw extruder (processing LIFE) have been 
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investigated It is to be noted that input data for and n 
variations are taken from Karwe and Jaluria (1990) and have been 
listed in Table-1 The input data for N variations are taken from 
Fenner (1979) and have been listed in Table-2. 


7.1 Effect of dimensionless flow rate (q^) 

Fig. 7 1 and Fig 7.2 show the effect of variation of i e , 
the flow rate (keeping all other variables such as frequency of 
screw rotation, screw geometry, etc. , constant) on the pressure 
gradient and pressure along the downchannel direction 
respectively. It is seen from Fig 7 1 that for q^ = 0 3, 0 35 and 
0.4, the pressure gradient is positive and is decreasing 
monotonically. For q^ = 0 45, the pressure gradient is slightly 
negative indicating that the drag flow situation has arrived. 
Also, as q^ increases, pressure gradient decreases. Consequently, 
from Fig. 7. 2 It is revealed that pressure rise decreases as the 
volumetric flow rate increases. The trend may be explained as 

follows. 


In the laboratory, various flow rates can be obtained by a 
variable resistance in the discharge end of the extruder. This 
accomplished by opening or closing the valve at the die end 
at low flow rates (l.e., d, = 0 3, 0,35 , 0 4, the pressure 
profiles increase nonotonlcally This is caused by the ability of 
the barrel to drag more poly^r than the net flow rate of the mil 
at each location In the extr^er channel. Therefore, the pressure 
gradient higher at lower flow rates. As flow rate is increased 



Pressure Gradient (3 Pqv/9z*) 


158 



Helical Length (Z^) 


Fig. 7.1 Pressure gradient distribution along the 
helical length for qy= 0 . 3 , 0.35,0.4, 0.45. 
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(q^«0 45 case) by opening a valve (at constant screw speed), the 
pressure gradient decreases, reaching a value at a certain point 
in the extruder where the drag flow generated by the advancing 
barrel equals the flow rate (net flow rate) No pressure gradient 
is develop)ed at this p>oint and further down the channel, negative 
pressure gradient is exhibited In the present case, since the 
metering section is the last section of the channel, no maxima is 
seen in the pressure curve and at z* = 0, the pressure gradient is 
already negative and remains negative till the channel exit 
Hence, q^ = 0 45, the pressure continues to decrease up to the die 
end (Fig 7.2). That is why, a value for q^ between about 0 4 and 
0.5 corresponds to the drag flow situation, where no pressure 
gradient is generated. For this reason, the limiting value of q^ 
on the higher side is taken as 0 4. 

Fig 7.3 shows the Nusselt number (NUjj) distribution along the 
downchannel direction for q^ = 0 3, 0 35, 0 4, and 0.45 As 

expected, Nusselt number is higher for higher flow rates. This is 
because, with rise in flow rate, the residence time distribution 
of the polymer material decreases which means that the fluid is 
heated less and therefore, the fluid temperature is lower. 
Therefore Nusselt number (Nu^) being indicative of the temperature 
gradient at the barrel surface increases. 

Fig7.4 depicts the bulk temperature distribution 

along the downchannel direction for q^ * 0 3, 0.35 , 0.4, and 0.45. 
As expected, the bulk temperature at any z»-location is higher for 
lower volumetric flow rates. The reason has been already stated in 
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Fig 7.3 Nussett number distribution along the helical 
length for qv=0.3, 0.35, 0.4, 0.45. 
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Fia 74 Bulk temperature distribution 
® Wical length for q, = 0.3, 0.35,0.4, 0.45, 
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table-1 


Parameter 

Value 

Power law index (n) 

0 5 

Screw helix angle (^) 

16 54° 

Griffith number (G) 

10 0 

Peclet nvimber (Pe) 

5000 0 

Temperature co-efficient of viscosity (p) 

1.61 


7.2 Effect of power-law index (n) 

Fig 7.5 and Fig 7.6 show that the effect of variation of 
power-law index n, on the pressure gradient OpVdz*) and pressure 
(p») along the downchannel direction respectively It is to be 
noted that n=0 5 and n=0.3 represent different grades of LDPE 
(i e. , of different melt flow index, MFI) The viscosity of the 
fluid having n=0.3 is greater than that for n=0.5. On the other 
hand n^l-0 represents a Newtonian fluid. It is evident that 
viscosity of a Newtonian fluid does not depend on the strain-rate. 
Furtherffiore. the viscosity of a Newtonian fluid (e.g. lubricating 
oil) is several orders of magnitude lower than that of molten 

polymers. 

inspection of Fig 7.5 reveals that pressure gradient for 
nai O is lowest in comparison with that for n*0.3 and 0.5. This is 
for n^l.O, viscosity is lowest and therefore, lowest 


becaus!® 
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Fig. 75 Pressure gradient distribution olong the 
helical length for n » 0.3, 0.5, 1.0 . 
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pressure gradient is generated For the aforementioned reason, it 
is seen from Fig 7 6 that pressure rise from inlet to outlet is 
naximum for the fluid with n=0 3, followed by fluids wi.th n=0 5 
and n=l 0 respectively 

Fig. 7. 7 depicts the Nusselt numbers, NUjj variations along the 
downstreeun direction of the channel. Nusselt number is lowest for 
the Newtonian fluid because of higher fluid temperature, the 
reason being greater strain rate leading to overall greater 
viscous heating although viscosity is low The strain rate is 
large because the fluid moves with greater velocity as the 
viscosity is smaller with respect to fluids with n=0 3 and n=0 5. 
The net effect is greater viscous heat dissipation (also evident 
from the term G 7 *^*^ in the non-dimensional energy equation of 
the fluid, where for Newtonian fluid n=1.0 as compared to n<l for 
other two fluids, G remaining same for all three). It is also 
noticed that NUjj becomes negative for the Newtonian fluid earlier 
in the downchannel length meaning greater heating of the fluid so 
that the fluid temperature near the barrel becomes greater than 

the barrel temjjerature. 

Fig. 7 . 8 shows that as expected, for the Newtonian fluid the 
bulk tenperature rise Is highest for the reasons sentloned In the 

previous paragraph. 


7.3 Effect of frequency of screw rotation m 


The frequency 


of screw rotation 


is the single most important 
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Fig. 7.7 Nusselt number distribution along the helical 
length for nsO.3, 0.5, 1.0 
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Fig 7.8 Bulk temperature distribution along the 
helical length for n = 0.3, 0.5, 1.0 
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variable, because It can be controlled best and is also most often 
varied The input data is listed in Table-2 Note that for this 
parametric study, data for a different extruder handling a 
different grade LDPE (i e. , with different power-law index, n) 
have been used Three screw speeds of 90, 100, 110 r.p.m have 
been selected because usually most plasticating screw extruders 
operate in this range The dimensional flow rate Q is kept 
constant. (Note that, in the laboratory this can be achieved by 
means of a valve at the die-end) However, in this case, will 
be changing since N is changing (meaning is different for 

different N’s) 


Fig. 7. 9 shows the pressure gradient variations along the 
downchannel direction for N = 90, 100, 110 r.p m . It is readily 
seen that lower screw speed results in lower pressure gradients 
and hence lower pressure rise in the extruder (Fig. 7. 10) and 
vice-versa This is because at higher screw speeds the barrel 
drags more amount of fluids but since constant Q is maintained the 
valve opening is reduced at the channel exit resulting in higher 
pressure gradient and hence, higher pressure rise 


Fig. 7. 11 shows Nusselt number (Nu^) variations along the 
downchannel direction for various screw speeds It is interesting 


to note that for N = 90 r.p.m , Nu^ at the channel exit is still 
positive unlike those for N = 100 and 110 r.p.m. in which cases 
NUjj’s are negative. This is because for lower r.p.m. (N = 90. in 
this case) viscous heat dissipation is less and hence the fluid 


temperature near the barrel 


at the die end has not reached a value 
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FiQ 7.9 Pressure gradient distribution along the 
helical length for Ns90, 100, 110 
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Fig. 7.11 


Nusselt number distribution along the 
helical length for N= 90, 100, 110 
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greater than or equal to the barrel temperature in sharp constrast 
to other two cases where viscous dissipations are higher 

Finally, Fig 7 12 shows the bulk temperature (0, ,, ) rise 

bulk 

along the extruder, which, as expected, indicates higher value for 
higher screw speeds, the reason being greater viscous heat 
dissipation for larger screw speeds 

TABLE 2 

Parameter Value 

Barrel inner diameter, D 

Screw channel depth, H 

Flight width, e 

Axial length, L 

Screw helix angle (<p) 

Flow rate, Q 

Power law index, n 

Reference viscosity, 

Reference strain rate, 

Temperature coefficient of viscosity, b 

T -T 

b(ref) 1 

T.-T 
1 0 

Specific heat, C 
Density, p 

Thermal conductivity, 

Thermal conductivity, 


120mm 
6mm 
12mm 
480mm 
17 66° 

510 kg/hr 

0 36 

10 8 kNs/m^ 

1 0 s"^ 
0.022 C“^ 

10 C 

10 C 

2000 J/kg K 
990 kg/m^ 

0 3 W/m K 
45 0 W/m K 



Bulk Temperature (Opulk^ 
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Fiq 7.12 Bulk temperature distribution along the 
helical length for N = 90, 100,110 
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CONCLUSIONS 



CHAPTER 8 


CONCLUSIONS 

A detailed numerical study of thermal transport processes In 
the metering section of a single-screw plastlcatlng extruder 
processing non-Newtonian power-law fluids has been carried out 
Unlike In the earlier works, conduction from the hot melt to the 
screw body has been taken Into account. In the process, three 
computer models of Increasing complexity have been developed and 
results based on each compared and contrasted with each other 

The main achievement of this work is the development of a 
quasi three-dimensional computer model that includes the Important 
cross-convection (1 e. , convection normal to the screw flights) 
thus leading to a very realistic prediction of screw surface 
temperature and pressure distribution in the downchannel 
direction. The nature of the screw surface temperature profile 
matches very well with the trend observed in experiments conducted 
by earlier researchers in this field. This clearly establishes the 
superiority of the present three-dimens Inal model in its ability 
to predict flow, thermal and pressure fields in the screw channel 
over earlier two-dimensional models. However, the two most 
important outcomes of the comparative study of the quasi-two 
dimensional and fully two-dimensional models also developed by the 
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present author are: (a) development of an efficient and novel 
numerical technique for computing the pressure graients accurately 
in the case of the fully two-dimensional model and (b) showing 
that the flow becomes locally fully developed at a short distance 
from the Inlet, - thus confirming the earlier use of fully 
developed velocity profiles in the quasi two-dimensional model as 
the inlet boundary condition to be not physically imrealistlc The 
main limitation of both the models is that they do not take into 
account the thermal convection normal to the screw flights and 
thus fall to predict realistically the screw surface temperature. 
The necessity of realistic prediction of screw surface temperature 
as well as pressure distribution is the motivation behind the 
development of the quasi-three dimensional model 

The velocity vector plots in the cross-sectional planes of 
the channel at various downstream locations reveal the circulation 
normal to the screw flights is much stronger as compared to that 
normal to the screw root. Consequently, thermal convection normal 
to the screw flights plays a major role in the mixing of the fluid 
and hence, in bringing the screw surface temperature very close to 
the barrel temperature This is also confirmed by the isotherm 
plots in the cross-sectional planes at various do^mstream 
locations which show almost uniform temperature near the exit of 
the channel indicating that the melt is throughly mixed. 

The parametric study based on the quasi three-dimensional 
model reveals that for values of dimensionless volumetric flow 
rate (q^) greater than 0.4, pressure gradient (dp»/dz*) is 
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negative (close to zero) which is not desirable because the very 
purpose of raising the fluid pressure from inlet to outlet will be 
defeated Therefore, the limiting value of should be 0 4 on the 
higher side It is also seen that for Newtonian fluid (n=1.0) 
pressure rise is low and bulk temperature rise is high as compared 
to two other grades of LDPE (n=0 3 and n=0.5). Also, higher screw 
speed results in higher pressure rise and higher bulk temperature 
rise along the downchannel direction of the extruder and 
vice-versa. 

The present quasi three-dimensional model has also formed a 
firm basis for a realistic design and optimization study of 
single-screw extruders 
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NON-UNIFORM GRID IN THE X-DIRECTION 



APPENDIX - A 


Non-uniform control volume face spacing is used in the 
x-dlrection The spacings are symmetric with respect to the 
mid-width location of the channel width Therefore, while 
calculating W/2 is used (see Fig Al). 
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Fig A1 Non uniform grids in X- direction. 


NI = Total number of grid points. 

NF = NI - 1 = Number of control volume faces 

NCAL * NF/2 

FI = Expansion factor (1 a. , >1) for gtW “‘'P 
x-direction 
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Non-uniform control volume face spacing is used in the 
x-direction The spacings are symmetric with respect to the 
mid-width location of the channel width. Therefore, while 
calculating Anjj, W/2 is used (see Fig.Al), 
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Fig.Al Nonuniform grids in X-direction. 


NI = Total number of grid points 

NF = NI - 1 = Number of control volume faces 

NCAL = NF/2 

FI = Expansion factor (i ©• . >1) for grid step length in 
x-direction 
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= Grid spacing adjacent to the screw flight 
As shown in the Fig 


NCAL-1 


= Ai) 


P NCAL-1 
1 ~ ^ 
1 FI - 1 


NCAL-1 


+ An, 



NCAL-1 
FI - 1 



NCAL-1 1 


2 


(Al) 


Let, A = 


, NCAL-1 , NCAL-1 

1 ■ ^ ^1 
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FI - 1 2 
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Therefore, 


AtJj = (W/2)/A 


(A3) 
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APPENDIX B 


HANDLING OF BOUNDARY CONDITIONS IN THE X-DIRECTION 


The temperature boundary conditions at x=0 and W/H are. 



To discretize it, a parabolic profile such as , 

2 

0=a+bx+cx 


(Bl) 


(B2) 


is assumed 

At X = x(l) = 0, (see Fig A2) 
0 = 0 , 


l.J 


Therefore, 


a = 0 


l.J 


At X = x{2) , 


0=0 , = e, , + bxC2) + cx(2) 

At X = x(3) , 

0=0=0, , + bx(3) + cx(3)^ 

j J A j J 

Solving Eq. (BS) and (B6), 

®2 a x(3)-^x(2) 

^ = {x(3) - xC 2 )}xl 2 )x(Jr ■ l.J x(2)xl3y 


Differentiating Eq (B2) w r t x. 


= b + 2cx 
dx 


At X = 0. 


de 

ax 


b 


(B3) 

(B4) 

(B5) 

(B6) 

(B7) 

(BS) 

(B9) 


184 




M\nd5o 


th 

^3C-2 

i>2eA 




Date Slip 


11785C 


This book Is to be retyrned on the 

date last stamped 




